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This  abstract  is  approved  by: 
Dissertation  Director: 


Julio  Gea-fianacloche 


This  dissertation  explores  the  connection  between  two  models  for  a  two-level  atom 
interacting  with  a  single  mode  electromagnetic  field.  The  first  is  a  semiclassical  model  in 
which  the  field  is  treated  classically  and  the  atom  is  treated  quantum  mechanically,  while 
the  second  is  a  full  quantum  model  in  which  both  the  atom  and  the  field  are  treated 
quantum  mechanically  (sometimes  known  in  quantum  optics  and  other  fields  as  the  spin- 
boson  model).  A  unique  pair  of  states,  which  remain  factorizable  for  long  times  and  are 
termed  quasiclassical  states,  are  derived  from  the  semiclassical  model.  The  quasiclassical 
states  are  found  to  match  the  evolution  of  the  quantum  model  well  for  fields  with  the 
average  photon  number  sufficiently  large.  These  quasiclassical  states  are  also  derived 
from  the  eigenvalues  and  eigenvectors  of  the  full  quantum  mechanical  model  alone. 

The  semiclassical  system  has  been  shown  to  exhibit  chaos.  We  identify  the  causes  of 
the  transition  to  chaotic  behavior  in  the  semiclassical  system  in  the  vicinity  of  the 
quasiclassical  orbits  and  examine  the  quantum  system  for  signatures  of  that  transition 
Several  authors  have  claimed  to  have  found  such  signatures,  but  by  using  the 
quasiclassical  states  we  show  these  supposed  signatures  are  not  unambiguously  related  the 
appearance  of  chaos  in  the  semiclassical  model.  The  quasiclassical  states  also  establish  a 
basis  for  comparing  results  for  calculations  which  use  the  rotating-wave  approximation 
(RWA)  with  results  from  those  which  do  not. 

Finally,  we  extend  the  work  by  considering  the  effects  of  up  to  three  atoms 
interacting  with  the  quantum  field.  The  result  is  that  the  system  behaves  more  classically 
than  for  a  single  atom  and  some  possible  signatures  of  chaos  are  identified  but  cannot  be 
established  unequivocally. 
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CHAPTER  1 


REVIEW  OF  PREVIOUS  WORK 

The  spin-boson  system  consists  of  a  two-level  system  coupled  to  a  harmonic 
oscillator.  This  system  has  been  studied  extensively  because  it  is  a  relatively  simple, 
closed  interaction  of  two  quantum  objects  and  has  many  different  physical 
manifestations.  In  this  work  we  shall  consider  a  two-level  atom  interacting  with  a 
quantized,  single  mode  electromagnetic  field, *>2  interacting  via  the  dipole 
approximation.  The  Hamiltonian  for  this  system  takes  the  form 

H  =  Ti(oa^a  +  \ticoQCj^  -\-tiga^{a+a^),  (1.1) 

where  (o  is  the  field  frequency,  coq  is  the  energy  spacing  of  the  two  atomic  levels,  g  is 
the  coupling  constant,  a  and  r/t  are  the  boson  annihilation  and  creation  operators,  and 
cr^  and  cr^  are  the  Pauli  spin  operators.  The  first  term  represents  the  energy  in  the  field, 
the  second  term  the  energy  of  the  atom,  and  the  third  term  is  the  interaction  energy. 
However,  other  relevant  physical  forms  include  a  two-level  system  such  as  an  atom, 
electron,  or  exciton  in  a  solid  coupled  to  a  phonon  mode  (see  for  example.  Refs.  3  and 
4);  an  ion  in  a  trap,  where  the  internal  excitation  of  the  ion  forms  a  two-level  system 
and  the  motion  of  the  ion’s  center  of  mass  is  described  by  the  oscillator;^’^  or  a  double 
well  system  under  the  influence  of  a  periodic  driving  force.^  There  are  more 
realizations  of  the  spin-boson  system,  but  these  few  provide  a  demonstration  of  the 
broad  application  of  the  problem. 

The  purpose  of  this  work  is  to  study  the  relationship  between  the  semiclassical 
and  fully  quantum  mechanical  versions  of  the  problem.  The  semiclassical  version,  in 
which  the  field  (oscillator)  is  treated  as  a  classical  system,  was  found  to  exhibit  chaos 
in  1976.*  Since  then,  the  spin-boson  system  has  been  studied  by  numerous  authors  in 
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an  attempt  to  discover  the  relationship  between  classical  chaos  and  “quantum  chaos.” 
In  what  follows,  the  remainder  of  Chapter  I  will  review  the  relevant  work  done  by 
others.  Chapter  2  will  examine  the  origins  of  chaos  in  the  semiclassical  model,  and 
Chapter  3  will  develop  and  examine  the  quantum  system  using  a  special  basis;  the 
quasiclassical  states.  Finally,  Chapter  4  will  use  the  quasiclassical  states  to  examine 
the  impact  of  semiclassical  chaos  on  the  fully  quantized  system. 

1.1  Semiclassical  Model — Prescribed  field 


Autler  and  Townes^  studied  the  semiclassical  version  of  Eq.  1.1  with  a 
prescribed  field.  Equation  1.1  becomes  for  this  case 

^ AT  =  2  cfy;  cos  oyl  (1.2) 

where  n  is  the  average  photon  number.  Because  this  is  a  Hamiltonian  periodic  in 
time,  Floquet’s  theorem  can  be  employed  to  find  two  solutions  of  the  form 

I  V/.('))  =  f  i;  A^e"‘'“\a)  -  ± 

\A:=-oo  l=-oo  J 

k.W>=-f 


where  \a)  and  |Z>)  represent  the  upper  and  lower  eigenstates  of  cr^,  respectively.  This 


is  the  reverse  of  the  notation  used  by  Autler  and  Townes,  but  has  been  used  to  provide 
consistency  with  other  publications’®'*^.  However,  the  A  and  B  coefficients  are 
unchanged  (except  for  the  trivial  relabelling  of  indices  k  and  /  from  including  only  even 
or  odd  integers  to  including  all  integers).  The  calculation  of  the  Floquet  characteristic 
exponent,  and  they4  and  B  coefficients  is  based  upon  infinite  continued  fractions 
and  given  in  detail  in  Ref  9.  The  two  special  solutions  (1.3),  which  shall  be  referred 
to  as  the  plus  and  minus  branches  henceforth,  will  appear  repeatedly  in  this  work. 

They  are  not  identified  specifically  in  Ref  9  but  come  from  Eq.  17  of  Ref  9  by  setting 
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Aq=0  and  ^0“^'  respectively.  (Incidentally,  there  is  a  sign  error  in  that  equation 
immediately  preceding  ^4.)  For  zero  detuning,  i.e.  <y=ri»o  ,  the  characteristic  exponent, 
and  the  A  and  B  coefficients  all  depend  on  a  single  parameter. 


This  parameter  is  the  ratio  of  the  Rabi  frequency,  2g-JjJ ,  to  the  nominal  resonant 
frequency,  coq,  and  characterizes  the  importance  of  the  counter-rotating  (energy  non¬ 
conserving)  terms  in  the  Hamiltonian  (1.1). 

Shirley  also  studied  the  semiclassical  problem  with  a  periodic  prescribed  field. 

He  approached  the  problem  using  Floquet  theory  instead  of  continued  fractions.  He 
expanded  the  periodic  Hamiltonian  using  complex  exponentials  and  he  solved  the 
problem  using  matrix.  The  approach  used  here  in  Chap.  2  will  be  similar;  however,  we 
will  go  beyond  his  work  by  considering  the  effect  of  the  atom  reacting  upon  the  field. 
Without  the  atomic  reaction,  chaos  is  impossible. 

1.2  Semiclassical  Model — Dynamic  field 

Although  numerous  authors  have  done  work  on  this  model,  I  will  review  in  this 
section  only  those  whose  efforts  have  directly  formed  a  basis  for  this  work.  The 
prescribed  field  model  cannot  admit  chaotic  solutions  because  it  is  exactly  solvable,  as 
seen  above.  Applying  a  different  approximation,  the  rotating  wave  approximation 
(RWA),  with  a  dynamic  field  cannot  admit  chaotic  solutions  either,  because  it,  too,  is 
exactly  solvable.  The  RWA  involves  dropping  the  two  terms  from  the  interaction 
portion  of  the  Hamiltonian  (1.1)  which  do  not  conserve  energy.  Replacing  cr^  with 
(a+  erf),  where  cris  the  spin  lowering  operator,  the  excluded  terms  are  aa  and  . 
These  terms  correspond  to  changing  the  atom  from  the  upper  to  the  lower  state  while 
annihilating  a  photon  and  raising  the  atom  from  the  lower  to  the  upper  state  while 


3 


creating  a  photon,  respectively.  Neglect  of  these  terms  generally  leads  to  only  very 
small  quantitative  differences.  However,  Tavis  and  Cummings'"’  pointed  out  that  for 
sufficiently  strong  fields  the  RWA  breaks  down. 

The  authors  of  Ref  8,  Belobrov,  Zaslavskii,  and  Tartakovskii  (BZT),  were  the 
first  to  identify  chaos  in  the  semiclassical  version  of  Eq.  1.1.  Expressed  in  terminology 
used  here,  they  found,  at  resonance,  the  entire  phase  space  was  chaotic  for  ^  n  132, 
or  0.35.  They  found  the  extent  of  the  chaotic  region  was  narrow  for  small  values 
of g’  and  increased  with  increasing  g.  The  chaotic  region  was  identified  by  examining 
the  time  evolution  of  the  trajectories  as  well  as  the  Euclidean  separation  between 
initially  close  trajectories.  They  attributed  the  chaos  to  the  perturbing  effects  of  the 
terms  excluded  in  the  RWA.  Later,  a  similar  result  was  published  by  Milonni, 
Ackerhalt,  and  Galbraith  (MAG)  on  a  slightly  different  model  using  Lyapunov 
exponents  and  Fourier  analysis  to  identify  chaos. However,  they  focused  on  a 
different  initial  condition  and  a  different  part  of  the  parameter  space  than  we  have. 

They  considered  an  initially  weak  field  with  a  large  number  of  interacting  atoms.  For 
weak  coupling  (small  g)  these  conditions  give  accurate  results  with  the  RWA.  In 
order  to  obtain  the  strong  fields  needed  to  cause  breakdown  of  the  RWA  given  the 
initially  weak  fields,  they  required  a  large  initial  population  inversion.  The  system  then 
had  enough  energy  to  generate  fields  strong  enough  to  cause  breakdown  of  the  RWA 
and  chaos. 

Fox  and  Eidson'^  compared  the  results  of  MAG  and  BZT  and  identified  the 
source  of  the  difference  in  the  models.  The  MAG  model  used  an  effective  polarization 
density  for  the  interaction  between  the  atom  and  the  field,  while  the  BZT  model  used 
an  effective  current  density.  Fox  and  Eidson  also  showed  the  effective  current  density 
arises  naturally  by  assuming  factorization  of  the  atom  and  field  in  the  fully  quantum 
mechanical  model.  Next,  they  demonstrated  by  a  change  of  variables  that  the  BZT 
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model  in  the  RWA  is  equivalent  to  a  classical  spherical  pendulum.  In  addition,  the 
initial  conditions  used  in  the  previous  studies  imply  the  pendulum  is  in  near-separatrix 
motion,  while  adding  the  non-RWA  terms  is  approximately  equivalent  to  providing  a 
periodic  perturbation.  Thus,  the  BZT  model  contains  a  periodically  perturbed,  near- 
separatrix  pendulum,  a  universally  accepted  chaotic  system.’"^  Moreover,  the 
numerical  results  showed  remarkable  similarity  between  the  pendulum  motion  and  the 
MAG  results,  despite  the  differences  in  those  systems. 

Kujawski  and  Munz’*^  published  a  study  of  the  frequencies  which  appear  in  the 
BZT  model.  Like  the  other  authors,  they  considered  initial  conditions  with  weak  fields 
and  near  full  inversion.  However,  they  found  the  frequencies  appearing  in  the  spectra 
were  combinations  of  harmonics  of  two  frequencies  which  they  identified.  They  noted 
the  values  of  these  frequencies  depended  upon  the  strength  of  the  field  in  the  cavity, 
but  they  did  not  specifically  identify  the  dependence. 

1.3  Fully  Quantum  Mechanical  System 

Ultimately,  most  of  the  interest  in  the  semiclassical  model  has  been  correlating 
the  chaos  observed  with  the  evolution  of  the  fully  quantized  system.  There  has  been  a 
great  flurry  of  activity  in  the  study  of  “quantum  chaos,”  despite  (or  perhaps  because 
of)  the  lack  of  a  generally  accepted  definition.  Because  Schrodinger’s  Equation  is  a 
first-order,  linear  equation,  chaos  in  the  strictest  sense  is  excluded.  However,  many 
classical  systems  have  quantum  mechanical  counterparts,  and  it  is  the  correspondence 
between  these  systems  that  forms  the  basis  of  “quantum  chaos.” 

At  least  two  distinct  approaches  have  been  used  to  investigate  this  area  for  the 

spin-boson  system.  First,  some  authors,  such  as  Grigolini,  have  studied 

dynamic  variables  such  as  the  purity  or  entropy  and  quantum  uncertainty  of  the 

system.  Their  results  were  based  on  two  general  initial  conditions,  the  atom  in  the 
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upper  state  and  the  field  (i)  in  a  coherent  state  with  a  small  average  photon  number  of 
«=5-10  or  (ii)  in  an  eigenstate  of  the  number  operator  with  /?=I0.  They  found  the 
evolution  of  the  quantum  system  was  marked  by  rapid  growth  in  the  field  uncertainty 
and  increase  in  the  entropy  of  the  atom-field  system  for  parameters  in  which  the 
semiclassical  system  was  chaotic.  They  attributed  these  results  to  the  effects  of 
quantizing  the  dipole.  Their  work  will  be  discussed  in  detail  in  Chap.  4. 

Another  approach  that  has  been  used  extensively  is  examination  of  the 
eigenvalues  of  the  full  quantum  Hamiltonian.  For  certain  types  of  systems,  such  as  a 
billiard  ball  in  a  stadium,  the  distribution  of  eigenvalues  correlates  closely  with  the 
chaotic  or  non-chaotic  behavior  of  the  classical  version  of  the  problem. For 
example,  classically,  a  ball  confined  to  a  circular  box  (stadium)  has  a  non-chaotic 
trajectory.  However,  if  the  stadium  becomes  an  oval  by  inserting  two  short,  straight 
parallel  line  segments,  the  trajectory  becomes  chaotic.  Quantum  mechanically,  the 
circular  stadium  has  many  eigenvalues  which  are  closely  spaced  together;  the 
distribution  of  nearest  neighbor  spacings  is  described  by  a  Poisson  distribution.  The 
eigenvalues  of  the  oval  stadium,  meanwhile,  exhibit  level  repulsion;  the  nearest 
neighbor  spacing  is  described  by  a  Wigner  distribution. 

Ku5  examined  the  spectrum  of  the  positive  parity  eigenvalues  of  Eq.  1  1  He 
found  the  eigenvalues  followed  neither  the  Poisson  nor  Wigner  distributions. 
Moreover,  the  eigenvalues  followed  a  very  regular  spacing,  even  for  large  values  of^ 
where  the  system  is  non-integrable.  Next,  he  measured  the  sensitivity  of  the  spectrum 
to  variation  in  parameters.  Sensitivity  to  variation  of  parameters  has  been  suggested 
as  a  characteristic  feature  of  quantum  chaos. He  found  only  the  lowest  lying 
eigenvalues  displayed  any  sensitivity — most  levels  were  quite  insensitive  to  changes  in 
(Oq.  He  concluded  from  the  spectrum  the  non-existence  of  “quantum  chaos,”  based 
upon  the  nature  of  the  two-level  system.  His  conclusion  was  supported  by  Graham 
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and  Hdhnerbach,25  who  reported,  furthermore,  that  the  nearest  neighbor  spacing 
distribution  does,  in  general,  approach  the  Wigner  distribution  for  values  of  angular 
momentum  larger  than y=l/2  and,  in  particular,  for y=9/2.  The  use  of  higher  angular 
momentum  corresponds  to  the  interaction  of  the  field  with  a  larger  number  of  atoms, 
so  a  more  classical  behavior  is  expected.  Earlier,  these  authors  presented  an 
approximate  expression  for  the  eigenvalues  developed  by  treating  the  o)q  term  in  Eq. 
1.1  as  a  perturbation. The  approximate  expression  gave  accurate  results  in  the 
region  of  most  interest  to  this  work,  large  n  and  g.  However,  when  they  applied  their 
results  to  discuss  the  relationship  between  the  behavior  of  the  quantum  system  and  the 
onset  of  semiclassical  chaos,^^  their  conclusions  are  questionable.  The  problem  with 
their  work  will  be  covered  in  Chap.  4. 

Although  substantially  more  material  could  be  included  in  a  thorough  review  of 
this  field,  such  a  review  is  beyond  the  stated  purpose  here.  Rather,  I  have  included 
only  the  material  which  serves  as  the  foundation  of  or  for  comparison  with  the  main 
body  of  this  work.  However,  the  next  section  presents  a  number  of  additional  articles 
which  were  examined  during  the  course  of  this  research. 

1 .4  Additional  Work  Reviewed 

This  section  will  very  briefly  describe  some  of  the  other  work  which  has  been 
done  on  relating  the  semiclassical  and  quantum  versions  of  the  spin-boson  model. 
However,  the  papers  included  in  this  section  did  not  make  a  vital  contribution  to  the 
present  effort.  They  are  included  here  for  the  sake  of  completeness. 

Cibils  and  coworkers  published  a  series  of  articles  comparing  the  results  of  the 

RWA  and  non-RWA  models  in  the  quantum  and  semiclassical  domains.2^-^o  [p, 

28  they  obtain  classical  Hamiltonians  which  yield  the  semiclassical  equations  of  motion 

for  both  the  BZT*  and  MAG'^  systems.  They  then  demonstrated  that  the  MAG 
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system  could  not  be  obtained  from  a  quantum  Hamiltonian  with  spin  1/2,  but  could  for 
any  larger  value  of  spin.  Next  in  Ref  29  they  use  a  phase  approximation  to  study  the 
system  semiclassically  and  quantum  mechanically.  They  use  the  phase  approximation 
to  show  a  back  reaction  of  the  atom  upon  the  field  is  always  present  in  a  consistent 
treatment  of  the  expansion  in  1//7.  Then  they  show  that  the  eigenvalues  using  the 
phase  approximation  and  the  full  quantum  mechanical  treatment  yield  similar  level 
statistics  (nearest  neighbor  distribution  and  Dyson-Mehta  A3  statistic)  in  the  non- 
RWA,  but  not  in  the  RWA.  Reference  30  reports  on  their  studies  using  the  Husimi 
function.  By  scaling  /j=0.01  they  pushed  the  quantum  system  into  the  semiclassical 
domain.  Then  they  examined  the  zeroes  of  the  Husimi  function,  which  form  a  type  of 
Poincare  map  (see  Sec.  2. 1  for  an  explanation  of  Poincare  sections)  for  the  quantum 
system.  They  found  the  zeros  of  the  Husimi  function  were  regularly  distributed  when 
making  the  RWA  and  were  somewhat  scattered  when  not  making  the  RWA. 

Studies  have  also  been  done  on  the  form  of  the  ground  state. The  more 
recent  of  these  references  found  good  agreement  between  the  exact  ground  state  and  a 
variational  ground  state  based  on  the  sum  of  two  coherent  states  for  intermediate 
values  of  the  coupling,  g,  and  the  detuning,  (y-a)Q.  These  authors  continued  their 
investigation  by  examining  Poincare  sections  of  the  field  for  both  the  semiclassical  and 
quantum  cases. In  the  semiclassical  case  they  used  a  procedure  similar  to  the  one 
used  here  in  Chap.  2.  For  the  quantum  case,  they  turn  to  the  Husimi  density.  They 
identify  certain  features  shared  by  the  semiclassical  Poincare  section  and  quantum 
Husimi  density.  In  particular,  they  find  that  islands  of  stability  in  the  Poincare  section 
are  visible  in  the  Husimi  density. 

Steeb,  e/  al.  consider  the  closely  related  one-fermion,  one-boson  Hamiltonian.-^"* 
When  making  the  RWA  approximation,  they  find  the  nearest  neighbor  spacings  of  the 
eigenvalues  are  distributed  almost  uniformly  up  to  a  maximum  value.  This  contrasts 
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(as  it  should)  with  behavior  associated  with  “quantum  chaos,”  for  which  the 
distribution  of  nearest  neighbor  spacings  typically  goes  to  zero.  For  the  non-RWA  the 
spacings  form  two  groups,  but  still  the  distribution  does  not  go  to  zero  for  zero 
spacing. 

Following  the  work  done  in  Refs.  23,  26,  and  34,  Schmutz  considers  an 
approximation  to  the  eigenvalues  of  Eq.  1. 1  in  the  high  n  limit.^^  Using  the  unitary 
transformation  of  Ref  3,  he  reduces  the  spin-boson  eigenvalue  problem  to  a  purely 
bosonic  system.  Then  by  extending  the  Hilbert  space  of  the  boson  to  negative  ii  states 
he  develops  an  asymptotic  scaling  law  for  the  eigenvalues.  This  scaling  law 
analytically  estimates  the  eigenvalues  for  large  n  for  any  value  of  the  detuning  and 
coupling. 

Finally,  Refs.  7  and  36  discuss  localization  in  a  double  well  system  with 
tunneling.  They  found  it  was  possible  to  localize  a  particle  in  such  a  system  by 
application  of  a  periodic  driving  force.  We  were  able  to  observe  this  localization  in  the 
semiclassical  model  in  the  chaotic  portion  of  the  parameter  space. 
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CHAPTER  2 


ORIGINS  OF  CHAOS 

I  have  used  several  tools  in  exploring  the  presence  of  chaos  in  the  semiclassical 
system.  The  first  part  of  this  chapter  is  devoted  to  describing  the  implementation  of 
those  tools.  The  second  part  of  the  chapter  will  apply  those  tools  to  explore  the 
semiclassical  model  in  the  neighborhood  of  the  periodic  orbits  obtained  from  Ref  9. 

2.1  Tests  of  Chaos 

The  widely  accepted  definition  of  chaos  is  a  “sensitive  dependence  on  initial 
conditions.”  Sensitive  dependence  is  demonstrated  by  an  exponential  separation  of 
initially  close  phase-space  trajectories.  The  most  straightforward  way  to  look  for  this 
behavior  is  by  simply  calculating  the  Euclidean  distance  between  the  points  in  phase 
space  and  plotting  them  on  a  semi-logarithmic  scale.  The  Euclidean  distance  is  defined 
to  be 

^  =  (2.1) 

where  x  and  x'  are  the  phase  space  coordinates  and  n  is  the  dimension  of  the  phase 
space.  For  example,  consider  the  Henon-Heiles  system,  which  was  thoroughly 
analyzed  by  Benettin,  Galgani,  and  Strelcyn.^^  This  is  a  classical  chaotic  system, 
containing  mixed  chaotic  and  non-chaotic  regions  for  appropriately  chosen  parameters. 
Figure  1  shows  the  separation  of  two  sets  of  trajectories,  one  from  a  regular  region  (a) 
and  one  from  a  chaotic  region  (b).  The  chaotic  trajectories  separate  exponentially 
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until  the  separation  is  approximately  equal  to  unit  distance,  which  is  the  approximate 
size  of  the  available  phase  space  for  the  chosen  set  of  parameters  (energy  E=0. 1 3,  q2=- 
0.1,  P2=0.  Refer  to  Ref  37  for  definitions  and  details  on  the  Henon-Heiles  system). 
Separation  of  the  regular  trajectories,  however,  is  much  slower  and  clearly  not 
exponential.  (E=0. 13,  q2=  0.0,  P2=0.0).  In  Sec.  2.2  we  will  see  similar  behavior  in  the 
semiclassical  spin-boson  model. 

Closely  related  to  the  exponential  growth  in  the  Euclidean  distance  between 
trajectories  is  the  concept  of  the  Lyapunov  exponent.  A,  which  gives  the  rate  of 
exponential  divergence  of  the  trajectories  in  phase  space  Because  trajectories  which 
are  exponentially  separating  shortly  reach  the  bounds  of  the  available  phase  space  for  a 
space  with  finite  volume,  accurate  calculation  of  A  requires  periodically  resetting  the 
trajectories  to  a  small  separation.  One  can  then  accurately  determine  the  average 
exponential  separation  rate.  Reference  37  gives  a  detailed  explanation  of  resetting  the 
trajectories  and  how  to  determine  A.  Most  importantly,  the  process  converges  on  A 
only  after  many  iterations.  If  A=0,  the  sequence  of  calculated  values  decrease  linearly 
on  a  log-log  plot,  i.e.  they  follow  a  power  law  decay.  Again  returning  to  the  Henon- 
Heiles  system.  Fig.  2  compares  the  calculation  of  A  in  both  ordered  and  chaotic 
regions,  using  the  same  initial  conditions  as  Fig.  I . 

Another  common  indication  of  chaos  is  in  the  frequency  spectra  of  the  dynamical 
variables  in  the  problem.  Regular  spectra  are  characterized  by  a  few,  well  separated 
peaks.  Chaotic  systems  have  noisy,  broadband  spectra.  Using  the  same  initial 
conditions,  consider  the  fast-Fourier  transform  (FFT)  of  one  coordinate,  shown  in  Fig. 

3.  The  first  graph  shows  the  spectrum  for  the  regular  trajectory.  It  shows  a  few  fairly 
narrow  peaks.  On  a  linear  plot,  only  half  a  dozen  or  so  would  be  visible.  The  second, 
in  marked  contrast,  shows  a  virtually  continuous  spectrum.  Even  on  a  linear  graph, 
much  of  the  structure  would  be  visible.  It  is  impossible  to  identify  particular  peaks  as 
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in  the  first  graph.  In  addition  to  examining  the  character  of  the  spectrum,  it  can  be 
useful  to  study  the  changes  in  the  spectrum  with  changes  in  parameters  or  initial 
conditions.  In  the  next  section,  we  will  see  how  numerical  investigation  of  the 
spectrum  led  to  an  analytical  result  which  helps  to  understand  the  origin  of  chaos  in 
the  vicinity  of  the  periodic  trajectories  described  by  Eq.  1.3. 

A  third  technique  to  identify  chaos  is  by  examining  the  Poincare  section.  The 
Poincare  section  is  a  stroboscopic  plot  of  the  trajectory.  At  periodic  time  intervals  or 
when  one  of  the  coordinates  takes  on  a  particular  value,  two  other  coordinates  (or 
other  phase  space  variables)  are  plotted  against  each  other.  Reference  37  gives  an 
example  similar  to  Fig.  4a.  In  this  figure,  whenever  the  first  coordinate,  is  zero  and 
the  corresponding  momentum,  /?],  is  greater  than  zero  a  point  is  plotted  on  the  graph. 
The  same  two  trajectories  were  used  in  this  figure,  and  form  two  distinct  sets  of 
points.  One  set  is  scattered  throughout  the  phase  space  and  results  from  the  single 
chaotic  trajectory,  while  the  other  set  forms  an  oval  and  results  from  the  non-chaotic 
trajectory. 

The  choice  of  variables  for  the  Poincare  section  is  not  unique;  although  some 
choices  are  less  suitable  than  others.  Figure  4b  shows  the  result  for  a  different  section. 
Notice  this  section  is  not  filled  in  as  much,  having  less  than  400  points  for  the  same 
trajectory  that  produced  over  1 200  points  in  the  oval  in  Fig.  4a.  Thus,  I  used  trial  and 
error  to  develop  the  Poincare  sections  used  in  Sec.  2.2,  and  found  the  results 
presented  there  to  be  satisfactory. 

Three  distinct  tools  for  examining  and  quantifying  chaos  in  a  classical  system 

have  been  presented:  Euclidean  separation  of  trajectories  and  the  closely  related 

Lyapunov  exponent,  spectral  analysis,  and  Poincare  sections.  Each  of  these  will  be 

used  to  first  identify  and  then  study  chaos  in  the  semiclassical  version  of  the  spin- 

boson  model.  First,  though,  I  want  to  develop  the  equations  of  motion  for  the 

semiclassical  model  and  relate  these  to  the  results  of  Autler  and  Townes.^ 
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2.2  Development  of  Semiclassical  Model 


Before  endeavoring  to  develop  the  equations  of  motion,  I  will  clarify  the 
notation  being  used.  The  states  of  the  two-level  atom  will  be  denoted  as  \a)  and  |Z»), 

|a)  being  the  upper  state.  The  energy  difference  between  |a)  and  \h)  is  Ticoq.  The 

electromagnetic  field  has  frequency  co  and  except  as  otherwise  noted  will  start  in  a 
coherent  state  |a),  with  a  = 

As  mentioned  in  Sec.  1.2,  Ref  16  provides  a  clear  development  of  the 
semiclassical  equations  from  the  purely  quantum  system.  Following  the  work  done 
there,  we  begin  with  the  Hamiltonian  in  Eq.  1.1  and  Heisenberg’s  operator  equations 


of  motion: 

(2-2a) 

=  (2-2b) 

cr^  =^[<T^,//]  =  2^o-„(a  +  nr^),  (2.2c) 

d+a^  =  — [a+«^,//l  =  -/Vy(nf-a^),  (2.2d) 

iti  '• 

d-d^  =  —\a-a^,H^  =  -ico{a+a^)-2iga^.  (2.2e) 

iTi  '■  ■' 

where  the  07  are  the  Pauli  spin  matrices.  Now  we  make  the  definitions 

xs(c7,),  (2.3a) 

y^{cyy).  (2.3b) 

z^{c7,),  (2.3c) 


17 


(2.3d) 


a,  = 


{a  +  a^^ 
2  ’ 


{a-a<) 


2/ 


(2.3e) 


In  order  to  apply  these  definitions,  we  must  take  the  expectation  values  of  Eqs. 
2.2.  However,  terms  which  contain  products  of  operators  must  be  treated 
approximately.  In  general,  (za)  (-)(<^)-  This  relationship  is  an  equality  only  for 

those  states  |  ^)  =  |  y/)0|<I>),  where  \  y/)  is  the  atomic  part  of  the  wavefunction  and  |<t>) 


is  the  field  part.  Thus  we  are  motivated  to  find  “quasiclassical  states,”  quantum  states 
which  remain  factorizable  for  a  long  period  of  time.  Using  work  done  with  the  RWA 
as  a  basis,*®" '2  we  expect  to  find  such  states  for  a  single  or  a  few  atoms  and  average 
photon  numbers  in  the  realm  of  10-100.  Further  details  on  finding  the  quasiclassical 
states  will  be  deferred  to  Chapter  3;  for  now  let  us  assume  we  can  treat  the  states  as 
factorizable  and  see  what  results  can  be  obtained. 

Upon  replacing  operators  with  their  expectation  values,  Eqs.  2.2  become-^^ 


x  =  -«oT 

(2.4a) 

y  =  (o^x-Agza, 

(2.4b) 

z  =  Agya^ 

(2.4c) 

a,  =  mar. 

(2.4d) 

a^=-coa^ -gx 

(2.4e) 

These  five  equations  are  the  semiclassical  equations  of  motion  for  the  atom-field 
system.  Scaling  the  system  by  coQi-^t,  (oI(Oq-^h,  A{glcO(^a^-^a^,  and  4(g/<yo)<^2~^^2  ’ 
and  defining  /(<y<yo)  =  \6^pil'n,  we  obtain  the  system  implemented  for 

numerical  computations; 


(2.5a) 

1 

X 

1! 

(2.5b) 

2  =  ya^, 

(2.5c) 

18 


(2.5d) 


a,  =//nr,, 

<7,  =  -/iOf,  -—X.  (2.5e) 

This  system  is  the  same  as  used  in  Ref.  1 8,  except  for  the  sign  of  the  field,  which  was 
chosen  to  match  Refs.  10-12.  Having  derived  this  for  the  general  case,  in  what  follows 
we  shall  only  consider  the  resonant  case,  //=1,  unless  otherwise  noted.  In  these 
systems,  the  three  components  x,  y,  and  z  form  the  atomic  Bloch  vector, '  which  has 

a  constant  unit  length.  In  addition,  energy  is  conserved,  which  has  the  form 


=  — +/7(af  +al  +  l)  +— 
rico„  2  ^  ’  co^ 


referring  to  (2.4)  or 


E  =  fi{z  +  xa^)  +  ^(af +a;) 


(2.6) 


(2.7) 


referring  to  (2.5). 

It  is  straightforward  to  generalize  this  to  interactions  with  larger  numbers  of 
atoms  in  order  to  obtain  a  more  classical  system.  First,  replace  the  cr,  with  the 
corresponding  angular  momentum  operators  Jj  in  the  Hamiltonian  (1.1)  and  in  Eq.  2.2. 
Let  the  angular  momentum  quantum  number  j  -  N/2,  where  is  the  number  of  atoms 
interacting.  Next,  change  the  definition  of  x,  y,  and  z  to  be  the  expectation  of  the 
corresponding  angular  momentum  operation  divided  by  j.  Equations  2.4a-d  are  then 
unchanged  and  Eq.  2.4e  is 

a^  =  -QXi^-Agxj .  (2.8) 

By  scaling  g—)-g/ ,  o,  — > «, yj2 j,  and  a,  a, Eqs.  2.4  are  restored.  Thus, 

in  the  semiclassical  version  of  the  problem,  the  equations  are  unchanged  by  suitable 
scaling  in  going  to  a  larger  number  of  atoms.  Interestingly,  since  Jf  -  >  «2y 

and  £  (=  g^^^ /2ft))  remains  unchanged.  Quantum  mechanically,  of  course,  there 
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remains  a  tremendous  difference.  For  clarity,  I  will  defer  further  discussion  of  values 
of  j  larger  than  1/2  until  Chap.  4. 

In  addition  to  being  useful  for  numerical  calculations,  Eq.  2.5  is  conceptually 
helpful.  Noting  that  scales  like  1//7  for  fixed  s,  one  can  easily  see  how  the  reaction 
of  the  atom  on  the  field  is  diminished  with  increasing  77.  Furthermore,  the  scaling 
makes  the  magnitude  of  the  field  dependent  only  upon  e.  In  Eq.  2.4,  a,  =  at  t=0. 
This  becomes,  in  Eq.  2.5,  a,  =  Se.  Increasing  77  with  s  fixed  makes  the  system  behave 
more  like  the  system  with  the  prescribed  field,  since  the  magnitude  of  the  field  does 
not  change  while  the  reaction  decreases  as  1/77.  Therefore  comparison  to  the  work  in 
Ref  9  becomes  particularly  straightforward,  since  the  system  (2.5)  reduces  to  that  one 
for  fixed  e  and  77 -^oo.  As  mentioned  earlier,  we  expect  to  require  photon  numbers 
somewhat  greater  than  77=10  in  order  to  have  states  which  remain  factorizable  and 
provide  a  reasonable  comparison  between  the  quantum  and  semiclassical  models. 

The  evolution  described  by  Eq.  1.3  is  clearly  non-chaotic.  Each  branch  evolves 
quasiperiodically.  How  does  this  quasiperiodic  behavior  in  the  limit  77  — >oo  become 
chaotic  as  77  decreases?  The  simplest  manner  to  compare  the  evolution  of  Eq.  1.3 
with  the  numerical  integration  of  Eq.  2.5  is  to  express  the  expectation  values  of  cr, 
when  taken  over  the  branches  1 1//^>  and  |  y/>.  In  order  to  arrive  at  a  general 
expression,  consider  a  linear  superposition  of  the  plus  and  minus  branches 
\il^)  -  C^\tf/ +  C_\y/ _) .  Since  this  state  should  be  normalized,  it  is  convenient  to 

write  C+  and  C.  as  cos{yH)  and  respectively.  Then  taking  expectation  values 

gives: 


00  oo 

cos((2a?  +  l)rty) 

/7=0  k=--oo 

oO  OO 

-  sin  y  y^X^2k^2n'~2k  ^2k  +  ]^2n’-2k-\)  codin' clj(  -  2Aj) 

/?'=-oo  JS:=-oo 


(2.9a) 
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(2.9b) 


y  =  COSrY,  +  cosln  ((2//  +  l)<<v) 

rt=0  k  =  -oo 

-sin  r  Z  +^2,..An'  2,  Ji2fj'6)-2Ajsin{2n'(Vt-2Aj) 

n'=-oo  k  =  -co 


=  cosy 


00  OQ  W  1 

E(4  - 4.,)+2E  Z(44..!. -b„..b,,.,.,  „)cos2n  cot 


\k=-ao 

CO  CO 


(2.9c) 


-sinrZ  CO  s((2«'  +  l)6y-2/l„/) 


m'=:-0O  /f  =  -cO 


When  y=0  or  ;rthe  trajectories  are  periodic,  corresponding  to  the  plus  and  minus, 
respectively  branches.  Since  ^J(0  is  generally  not  a  rational  number,  almost  all  other 
trajectories  are  quasiperiodic.  Using  the  methods  of  Ref.  9,  we  can  plot  21^  vs.  £,  as 
in  Fig.  5  below. 

Autler  and  Townes  mention  the  problems  associated  with  multiple  solutions  for 
the  characteristic  exponent,  and  we  emphasize  it  as  well.  Examining  Eq.  1.3  one  can 
readily  observe  that  if  is  a  Floquet  exponent,  2n(Oi-A^  is  also  a  Floquet  exponent. 
This  requires  simply  factoring  out  e'^’^'^  from  within  the  summation  while  relabelling 
the  indices  k  and  /.  This  process  can  be  repeated,  giving  an  infinite  number  of 
solutions,  none  of  which  are  physically  distinct.  However,  replacing  with  -Ag  results 
in  swapping  the  plus  and  minus  branches.  (This  result  does  not  come  directly  from  Eq. 
1 .3  because  the  ^’s  and  B's  depend  on  Ag;  the  existence  of  a  second  set  of  solutions  is 
discussed  in  Ref  9.)  The  new  solution  is  now  physically  distinct.  Alternatively,  the 
two  branches  can  be  considered  as  representing  two  solutions  with  positive  and 
negative  Floquet  exponents. 

Before  proceeding  further,  consider  the  appearance  of  these  solutions  in  terms  of 
the  Bloch  sphere.  Figure  6  shows  a  sequence  of  plus  branch  trajectories  for  different 
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values  of  e.  The  number  printed  by  each  sphere  displays  the  value  of  e  for  that 
trajectory.  For  small  values  of  irthe  trajectories  remain  in  near-equatorial  orbits — the 
atom  experiences  little  inversion.  However,  increasing  £•  increases  the  inversion,  so 
that  for  £=0.368  the  trajectory  actually  passes  through  the  upper  atomic  state.  For  the 
choice  of  field  phase  given  in  Eq.  1 .2,  the  initial  point  in  all  of  these  plots  is  in  the  x-z 
plane,  with  x>0,  as  shown  by  the  small  dot.  For  £  >  0.5 1 5,  however,  the  initial  point 
has  X  <  0,  and  for  even  larger  e  it  switches  back  again,  but  always  with  z  <  0.  The 
minus  branch  trajectories  are  obtained  by  reflection  through  the  origin,  but  are 
otherwise  identical. 

The  evolution  of  the  full  semiclassical  system  given  by  Eq.  2.5  yields  very  similar 
trajectories  for  sufficiently  large  77,  using  for  the  initial  condition  the  point  given  by 
Eq.  1.3  at  f=0.  It  is  convenient  to  define  the  starting  point  in  terms  of  the  spherical 
angles  in  the  Bloch  sphere.  In  that  case,  we  can  write 


I  W. (0))  =  +  sin  ^\b) 

\  W- (0))  =  sin  ^\a)  -  cos-^|A) 


(2.10a) 

(2.10b) 


where  6^  is  the  polar  angle  (colattitude)  measured  from  the  positive  z  axis. 
Comparison  with  Eq.  1.3  yields 


^  Jt  =  -00 

A  plot  of  0^  versus  e  is  shown  in  Fig.  7.  As  the  initial  point  passes  through  the  lower 
atomic  state,  becomes  greater  than  n.  The  initial  point  passes  through  the  lower 
atomic  state  twice,  once  at  £  =  0.5 1 5  and  again  at  £  =  0. 83 1 .  The  orthogonal  state  is 
located  at  the  opposite  end  of  the  diameter,  and  so  it  passes  through  the  upper  atomic 
state  at  the  same  values. 


25 


Differences  between  Autler-Townes  trajectories  and  the  semiclassical  evolution 
do  exist,  although  they  can  be  made  arbitrarily  small  given  sufficiently  large  Ti .  The 
first  effect,  noticeable  even  for  small  £,  is  an  interaction-induced  detuning.  Rather  than 
oscillating  at  the  frequency  coq,  the  trajectory  oscillates  at  <y',  which  is  slightly 

different.  By  writing  j£’(/)as 

00 

x(/)  =  ^je„cos(2w  +  l)<y7  (2.12) 

rt=0 

and  the  field  as 

o,(/)  =  a,(0)cos<y7  (2.13) 

(with  ai(0)=8£)then  inserting  into  Eq  2.5d,e  the  interaction  induced  detuning  becomes 

=  (2.14) 

2Vw  « 

by  matching  terms  with  the  same  time  dependence.  One  can  see  clearly  from  this 
relationship  that  the  detuning  vanishes  for  large  /7,  as  it  should.  Equation  2. 14  also 
shows  the  detuning  will  be  much  smaller  than  co  since  £  ^  1  and  Xq  <  1,  while  /7  >  10, 
Of  course,  the  precise  value  of  Xq  depends  on  the  detuning  used  to  calculate  it.  To 
obtain  even  better  agreement  between  the  Autler-Townes  and  semiclassical  models, 
this  procedure  can  be  iterated — finding  Xq  assuming  no  detuning,  then  using  the 
detuning  given  by  that  Xq  to  calculate  a  new  Xq.  Because  the  detuning  is  small,  this 
process  converges  after  no  more  than  two  iterations  for  f  ;S  1 . 

The  detuning  found  in  this  manner  gives  excellent  agreement  with  the  detuning 
found  by  numerically  integrating  Eq.  2.5  The  detuning  of  the  numerical  results  can  be 
found  most  easily  by  converting  the  field  expectation  values  into  a  complex  phasor 
a-a,  +  ia^ .  The  phasor  rotates  approximately  as  e"'®  '.  Multiplying  a  by  e'®'  cancels 
all  but  the  detuning,  which  can  then  be  determined  by  fitting  the  remaining  phase  drift 
to  a  straight  line. 
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The  detunings  for  the  two  branches  are  approximately  equal  in  magnitude  and 
opposite  in  sign.  This  is  obvious  from  Eq.  2.9a,  by  setting  /=  tt to  obtain  the  minus 
branch.  Clearly  Xq  for  the  minus  branch  is  opposite  in  sign  and  equal  in  magnitude  to 
Xq  for  the  plus  branch  if  we  neglect  the  detuning  in  Eq.  2.9a.  Furthermore,  since  the 
detuning  is  small,  the  effect  it  will  have  on  Xg  when  included  will  be  small  as  well. 

2.3  Chaos  in  the  Semiclassical  Model 

We  have  seen  the  semiclassical  model  has  periodic  trajectories  which  are  near 
the  Autler-Townes  solutions.  In  this  section  we  will  see  these  trajectories  are  stable 
for  sufficiently  large  n,  but  become  unstable  as  77  decreases,  resulting  in  chaos.  For 
some  values  of  s,  the  plus  branch  requires  a  larger  value  of  77  to  remain  stable  than  the 
minus  branch  requires.  For  other  values  of  s,  the  minus  branch  is  the  first  to  become 
unstable  as  77  decreases.  The  signature  of  the  instability  in  the  trajectory  is  found  in 
the  spectra;  As  77  changes,  the  location  of  the  peaks  shift  from  the  values  predicted  by 
the  Autler-Townes  solutions,  the  peaks  merge,  and  complex,  chaotic  spectra  result. 
First  we  will  examine  the  spectra  resulting  from  numerically  integrating  Eq.  2.5.  Then 
we  will  attempt  to  understand  the  numerical  results  analytically.  By  linearizing  the 
system  around  the  Autler-Townes  solutions  (2.9),  we  can  develop  an  approximate 
expression  for  the  frequency  shifts  and  understand  why  one  branch  becomes  unstable 
before  the  other. 

Now  consider  trajectories  which  are  near,  but  not  precisely  on  the  plus  and 
minus  branches.  Figure  8  shows  the  FFT  of  the  field  for  77  =  1 024  and  s  =  1  with  y=- 
.  1  ;r.  In  addition  to  the  prominent  peak  at  <y=l  .00  (recall  time  has  been  scaled  by  (Oq  in 
all  numerical  calculations),  additional  peaks  at  ^<>=0.13,  1.87,  2.13  and  3.00  are  visible. 
These  peaks  approximately  correspond  to  2<y-2/l„,  2/l„,  4rt>-2T^,  and  3a),  respectively. 
Of  course,  in  the  Autler-Townes  model,  the  only  frequency  appearing  is  O).  The  lower 
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part  of  the  graph  shows  the  equivalent  result  for  x.  Although  the  spectrum  is  much 
more  complicated,  most  of  the  frequencies  correspond  to  those  just  listed,  and  which 
do  appear  in  the  Autler-Townes  model  (see  Eq.  2.9).  However,  the  pair  of  peaks 
bracketing  co=l  and  the  four  peaks  bracketing  <y=3  (two  above  and  two  below)  do  not 
arise  from  the  Autler-Townes  model.  These  represent  higher  harmonics  and  are 
another  effect  of  the  atom  reacting  on  the  field.  Table  1  lists  the  frequency  of  each 
peak  from  Figure  8b  and  gives  the  combination  of  <y'  and  needed  to  produce  it. 

If  we  consider  a>'  and  each  as  fundamental  frequencies,  a  parity  for  each 
harmonic  can  be  calculated  by  summing  the  factors  multiplying  these  fundamentals. 
Notice  the  parity  in  Table  1 — the  total  parity  is  odd.  Even  harmonics  of  <y'  are  linked 
with  odd  harmonics  of  2/1^  and  vice  versa.  A  spectrum  of  z(/)  shows  the  same  type  of 
behavior,  but  with  only  even  parity  components.  The  parity  relationship  is  a  natural 
consequence  of  the  structure  of  Eqs.  2.5  and  2.9.  The  prescribed  field  contains  only 
odd  parities,  so  the  prescribed-field  solution  (Eq.  2.9)  for  x  and  7  contains  odd  parities, 
while  z  has  even  parities.  When  the  atom  is  allowed  to  react  on  the  field  (Eq.  2.5),  the 
parity  structure  may  (or  may  not)  be  preserved.  As  long  as  no  new  frequency 
components  are  introduced  by  another  mechanism,  the  structure  will  be  preserved. 
However,  if  frequencies  of  the  “wrong”  parity  are  introduced  by  some  means  into  one 
variable,  the  system  (2.5)  will  allow  those  frequencies  to  spread  to  other  variables.  We 
will  see  shortly  this  is  the  route  to  chaos  from  the  neighborhood  of  the  Autler-Townes 
trajectories. 

The  obvious  question  arises:  How  can  frequencies  of  the  “wrong”  parity  get 

introduced?  They  are  the  result  of  a  shift  in  various  frequency  components  due  to  the 

reaction  of  the  atom  upon  the  field.  Using  the  harmonics  given  in  Table  1  we  obtain 

agreement  to  four  decimal  places  between  the  calculated  harmonic  frequencies  and  the 

location  of  peaks  from  the  FFT.  However,  to  get  that  agreement,  we  must  use  the 

fundamental  frequencies  <y'  and  from  the  FFT,  not  the  Autler-Townes  model.  The 
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Harmonic _ Parity 


0,1316 _ 2a)'-2X„ _ 2-1=1 

0.4739  7co'-SX„  7-4=3 


.7372  3co'-4A,  3-2=1 


1.00045 


1.2637  5a) -4^  5-2=3 


1.3427  8co'-lOA„  8-5=3 


1.86926 _ 24 _ 1_ 

2.1325 _ 4co'-2A„ _ 4-1=3 

2.4748  5(o'-8A„  5-4=1 


2.7381 _ a)'-4A„ _ l-2=-l 

3.0014  30)'  3 


3.2646  7a)'-4X,  7-2=5 


3.3436  6(o'-\0X,,  6-5=1 


Table  1.  List  of  frequencies  found  in  FFT  of  x(t)  from  numerical  integration  of 
_ Eq.  2.5,  corresponding  harmonics,  and  parity. _ 
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Autler-Townes  model  predicts  the  frequencies  should  be  (o'=\  .0Q041  co  (using  Eq. 

2. 14)  and  2/l^,=  l  .86360ry.  Including  the  detuning  has  no  effect  on  the  calculated  value 
of  2A^,.  Thus,  while  the  Autler-Townes  model  predicts  the  detuned  co'  frequency  to 
within  0.002%,  it  does  not  account  for  the  discrepancy  of  0.3%.  One  is  tempted 
to  attempt  a  procedure  similar  to  the  one  which  led  to  Eq.  2. 14  to  find  a  detuning  for 
2A^.  Unfortunately,  such  a  simple  procedure  does  not  work.  The  difference  is  in  the 
nature  o^2X^.  It  arises  as  an  eigenvalue  in  the  solution  of  a  Floquet  eigenvalue 
problem,  and  so  it  involves  all  the  frequency  components  which  appear  in  the  problem, 
not  just  one. 

Let  us  examine  the  behavior  of  these  frequencies  with  changing  n .  Figure  9a 
shows  the  variation  in  the  frequency  of  the  2a)'-2X^,  harmonic  as  n  is  varied  and  shows 
no  obvious  pattern.  However,  plotting  the  same  frequencies  versus  1/77,  as  in  Fig.  9b, 
shows  a  distinctly  linear  dependence  for  sufficiently  large  77.  The  frequencies  from  the 
minus  branch  require  a  brief  explanation.  In  terms  of  the  components  of  the  Bloch 
sphere,  all  frequencies  appear  inside  of  cosine  functions  (see  Eq.  2.9)  and  so  the  sign 
of  the  frequencies  is  unimportant.  However,  Eq.  1.3  shows  the  actual  frequencies  of 
the  plus  and  minus  branches  are  of  opposite  sign.  In  order  to  compare  the  two 
branches,  the  absolute  values  of  the  frequencies  have  been  used.  Furthermore, 
because  the  frequencies  come  from  taking  FFTs  of  real-valued  time  series,  only  the 
positive  frequencies  are  significant.  Later  we  will  see  that  the  actual  sign  of  the 
frequencies  is  important,  but  for  now  we  are  primarily  interested  in  the  magnitudes. 

The  fit  to  the  plus  branch  intersects  the  co  axis  at  0. 1365  while  the  fit  to  the 
minus  branch  intersects  at  <y  =  0. 1372.  Both  provide  good  agreement  with  the 
Autler-Townes  model  prediction  of  2-1 .8636=0. 1364.  The  difference  is  attributed  to 
the  uncertainty  in  the  locations  of  the  peaks  in  the  FFT  and  in  fitting  a  straight  line  to 
the  locations  of  the  peaks. 
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In  order  to  find  the  correction  to  then,  let  us  return  to  the  original 
eigenvalue  problem  and  search  for  a  linear  dependence  upon  Mn  .  First  we  must 
express  the  system  (2.4)  as  an  eigenvalue  problem.  Let  us  consider  first  the  Autler- 
Townes  model  with  no  atomic  reaction  on  the  field.  The  reaction  comes  from  the 
term  gx  in  Eq.  2.4e.  Deleting  that  term,  we  can  write  a,  =  4^ cosoM ,  and  along  with 

Eqs.  2.4a-c,  we  obtain  a  system  of  linear  equations  with  periodic  coefficients. 
According  to  Floquet’s  theorem,  there  are  three  solutions  (since  there  are  three 
dimensions)  of  the  form 

X  =  ,  (2.15) 

n=~oo 

where  v  is  one  of  three  distinct  Floquet  characteristic  exponents.  In  addition,  y  and  z 
can  be  written  in  the  same  form.  Substituting  these  forms  into  Eq.  2.4a-c,  we  obtain 
an  eigenvalue  problem  which  can  be  truncated  and  solved  numerically.  In  doing  so, 
we  obtain  three  values  of  v={0,  2A^,,  -2A^}.  The  solution  v^O  gives  the  periodic 
solution  of  Eq.  2.9  (;^0,;r).  The  other  two  solutions  must  be  combined  with 
appropriate  choice  of  signs  to  obtain  real  solutions  which  correspond  to  Eq.  2.9  with 
y=  ±7dl.  Although  in  the  Hilbert  space  of  the  atomic  system  (i.e.  in  terms  of  |r/)  and 
|Z>))  the  plus  and  minus  branches  represent  orthogonal  states,  in  the  Bloch  sphere  they 

differ  only  by  a  minus  sign,  so  they  are  not  independent  solutions  of  the  eigenvalue 
problem  derived  from  Eq.  2.4. 

Now  to  obtain  the  correction  to  we  must  treat  the  atomic  reaction  term  as  a 
perturbation.  To  do  so,  we  shall  include  the  two  field  equations  in  our  eigensystem. 
However,  the  system  is  then  no  longer  a  linear  system.  If  we  restrict  our  interest  to 
trajectories  near  one  of  the  branches,  we  can  write  x  =  x'*^  +  5x,  where  5x  will  be 
small,  and  so  on  for  z,  aj  and  Oi-  Then  substituting  into  Eq.  2.5  we  obtain 
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x^^  +  a  =  -y^^  -dy 

(2.16a) 

y^^  +Sy  =  x'‘^  +Sx-z''^&J, 

(2.16b) 

z'^^  +  a  =  +y‘^^&l^  +a^^Sy 

(2.16c) 

-Px'*^  -  Sa^  - pdx 

(2.16d) 

+  5a^-  +  6a^. 

(2.16e) 

In  Eq  2. 16a-c,  the  terms  which  involve  only  the  unperturbed  solutions  cancel, 
since  those  are  the  unperturbed  Autler-Townes  equations.  The  field  equations  must  be 
handled  more  carefully.  Writing  Eqs.  2. 16d-e  as  a  single  second  order  differential 
equation,  we  get 

+  &i^  +  +Sa^-  +  &:) .  (2. 1 7) 

From  here  we  substitute  -  Stcos/y/,  where  fj=(o’l(OQ.  Then  using  the  expression 
for  <y' given  by  Eq.  2.14,  Eq.  2.17  then  becomes 

&r,  +  &f,  =  -pSx  -  -  Xg  cos fji) .  (2. 1 8) 

The  term  in  parentheses  in  Eq.  2.18  is  a  small  correction  to  an  already  small  quantity, 
but  cannot  he  a  priori  ignored.  However,  numerical  integration  of  the  system  with 
and  without  that  term  showed  negligible  differences  for  ii  in  the  linear  range  of  Fig.  9. 
Therefore,  the  term  will  be  dropped  and  we  obtain  the  linearized  semiclassical  system; 


Sx  =  -5y  (2.19a) 

Sy  =  Sx  -  5z  -  z'*^  Sa^  (2.19b) 

&  =  af^y+y^‘^5a^  (2.19c) 

Sa,=-5a^-pdx  (2.19d) 

Sa,^Sa,  (2.19e) 


While  this  system  appears  almost  identical  to  the  original  semiclassical  model,  all 
the  coefficients  are  now  periodic  because  we  have  expanded  about  the  unperturbed 
trajectory.  Furthermore,  the  system  can  be  expanded  as  a  sum  of  complex 
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exponentials,  such  as  (Sc  =  e"''  ^  Sx^e'"'“  .  Then  the  equations  of  motion  can  be 


written  in  two  parts:  the  unperturbed  part  (which  does  not  involve  p  or  8d)  and  the 
perturbation  (which  does  involve  /?  or  5d).  We  obtain,  in  matrix  form,  with  /j=(o'I(Oq 


(Sa  > 

-inju 

1 

0 

0 

0  ^ 

(Sa  \ 

00 

OQ 

-1 

-in/j 

0 

0 

0 

z- 

W 

li 

+ 

Ik 

0 

0 

-in/j 

-1 

0 

«=-00 

rt=-oo 

0 

0 

1 

-inn 

-8£C0S/// 

Sy„ 

J 

< 

0 

0 

0  %ecosfjt  -injj,  J 

( 

0 

0  0 

0  0^ 

j/(  v'+«//)r 


z 


0 

0 

yAT 


0  0  0 
0  0  0  0 
0  0  0  0 
0  0  0  0 


Sy„ 

Sz. 


(2.20) 


«  y 


Replacing  cosjjt ,  and  with  their  expansions  in  complex  exponentials  allows 
Eq.  2.20  to  be  written  as  an  eigenvalue  problem,  nc=(A+B)x.  The  actual  forms  of  A 
and  B  are  too  large  (5xoox5xoo)  to  be  included  here,  but  truncated  forms  are  included 
in  Appendix  A.  The  first  matrix.  A,  is  the  unperturbed  part,  which  is  satisfied  by  the 
Autler-Townes  solution.  A  specific  value  of  «  must  be  selected  in  order  to  calculate 
the  detuning.  Using  £=1.0  and  77=1024,  the  unperturbed  eigenvalues  when  expanding 
about  the  plus  branch  are 


-1.1378 

-0.1374 

0.8631 

-1.0009 

-4.7- 10-^ 

1 

... -1.0005 

0 

1.0005--- 

-1 

4.7-10'^ 

1.0009 

-0.8631 

0.1374 

0.8631 

(2.21) 


The  eigenvalues  come  in  sets  of  five,  with  each  set  separated  in  frequency  by  //.  The 
eigenvector  for  the  eigenvalue  v=0  corresponds  to  the  periodic  trajectory  given  by  Eq. 
2.9  with  y=0,n,  as  discussed  above.  Likewise,  the  eigenvalues  v=±2Xg  (=±0. 1374) 
have  already  been  discussed.  Note,  however,  that  the  eigenvalue  obtained  by 
following  the  convention  of  Autler  and  Townes  (see  Fig.  1  of  Ref  9  and  the  discussion 
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which  immediately  follows  it)  yields  2/iy<y=l  .8636=2//-k  Again,  care  must  be  taken 
to  avoid  problems  from  having  multiple,  physically  equivalent  solutions.  The  field 
components  of  both  of  these  eigenvectors  are  all  zero. 


The  new  eigenvalues,  v=±{p-\),  are  simply  the  interaction-induced  detuning,  S. 
The  corresponding  eigenvectors  give  the  periodic  field  oscillations.  Like  the  solutions 
using  ±2Ag,  the  periodic  field  solutions  must  be  obtained  from  linear  combinations  of 
the  eigenvectors  in  order  to  obtain  real  solutions.  Since  the  eigenvectors  are 
particularly  simple  for  these  eigenvalues,  consider  how  the  real  solutions  are  formed. 
The  &f|  part  of  the  eigenvector  for  the  positive  detuning  eigenvalue,  v=(j^-\)  can  be 
expressed  as 


Sr  ' 


((n+lV-l)/ 


(2.22) 


where  „  is  the  Kronnecker  delta.  If  v=m/.iA,  then  Eq.  2.22  can  be  generalized  to 


5a\=—j=y^Smy  - 

'  42^  '  42 

Similarly,  the  Soi  part  of  the  eigenvector  has  an  iH2  in  the  same  frequency 
component.  In  the  same  manner,  the  negative  detuning  eigenvalue  yields 


(2.23) 


i( /!//+(/;»//+ l))l 


^  ^  ^  ,7 


^2  -  4^^^ 


(2.24) 


4^ 


Therefore,  the  real-valued  Sa^({)  is  formed  by  the  sum  of  the  eigenvectors,  while 
Sa2(t)  is  formed  by  the  difference. 

The  change  in  the  Floquet  exponent  due  to  the  perturbation  (i.e.,  the  reaction  of 
the  atom  on  the  field,  described  by  the  second  matrix,  B)  can  be  found  using  standard 
perturbation  techniques.  Because  of  the  parity  considerations  discussed  at  length 
earlier  in  this  section,  first  order  perturbation  theory  yields  no  change  in  2A„.  The 
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(2.25) 


expression  for  the  first  order  change  in  an  eigenvalue  is 

where  W  is  the  perturbation  and  |  ‘ire  the  unperturbed  eigenvectors.  The 


perturbation  matrix  elements  connect  the  field  and  atomic  portions  of  the  eigenvectors. 
However,  no  eigenvector  has  non-zero  components  in  both  the  field  and  the  atomic 
portions.  Therefore,  the  first  order  perturbation  is  zero  for  all  the  eigenvalues. 

Care  must  be  taken  in  using  second  order  perturbation  theory  because  the 
perturbation  is  non-Hermitian.  The  appropriate  expression  for  the  second  order 
correction  to  the  wth  eigenvalue  is 


(2.26) 


We  are  interested  in  the  perturbation  to  the  eigenvalues  v=^2X^  (which  are  equivalent 
to  the  eigenvalues  v=+0.1374).  Since  these  correspond  to  atomic  eigenvectors,  based 
on  the  discussion  in  the  preceding  paragraph  the  only  eigenvectors  which  will 
contribute  to  the  sum  are  the  field  eigenvectors.  Let  the  unperturbed  atomic 
eigenvector  correspond  to  the  eigenvalue  v(=-0. 1374).  Let  the  field 


eigenvalues  be  denoted  by  or„,±=w/v±^  with  the  corresponding  eigenvector  denoted  by 
Having  defined  these  quantities,  Eq.  2.26  becomes 


(2.27) 


where  the  ±  indicates  summation  over  m  first  with  the  upper  sign  and  then  with  the 
lower  sign.  Further  simplification  is  achieved  by  noting  that  the  only  non-zero  part  of 
which  survives  is  the  <501  portion,  while  only  the  5^2  portion  of 

survives.  The  non-zero  component  from  eigenvector  is  -(wil),  as  discussed 

above.  In  order  to  find  the  matrix  elements  of  B,  first  express  Eq.  2.9  in  terms  of 
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complex  exponentials.  Then  let  x„,  and  z„  be  the  coefficients  of  with  and 
and  2',,  be  the  coefficients  of  with  y=Kl2.  Then  the  matrix  elements  of 

B  are  given  by  x„,y„,  and  z„.  The  components  of  lyz+^O))  can  be  obtained  either  from 

the  eigenvector  corresponding  to  v  or  from  x'„,  y'„,  and  2',,  after  they  have  been 
normalized  to  unity.  Equation  2.27  can  be  expressed  as 


I(x‘ 


^n±\+ni  ^rt±l+wj 


v-m/j+S 


(2.28) 


Finally,  the  summation  in  the  numerator  of  Eq.  2.28  is  exactly  equal  to  -ix'^:^_„,.  Thus 
the  second  order  perturbation  can  be  written  in  the  fairly  simple  form 


m 


v-m^+S 


(2.29) 


Fortuitously,  even  for  values  as  large  as  £=1,  only  a  single  term  from  each  sum 
over  m  is  significant,  m=+\  from  the  summation  over  the  upper  signs  and  m=-\  from 
the  summation  over  the  lower  signs.  This  simplification  leads  to  the  final  result 


Av/^>  = 


2AM«  _ 

•^0  _ 

1 

1- 

2//  “  \  \ 

2ju 

^0  J 

1  "0  J 

(2.30) 


where  x'q  is  the  component  of  the  eigenvector  which  oscillates  at  the  frequency  2<y - 

2A^.  Thus  we  obtain  the  desired  \/n  dependence  for  the  variation  of  2 Forf=l, 

Eq.  2.30  evaluates  to  approximately  6.9//7 ,  This  value  is  in  reasonable  agreement  to  a 

linear  least  squares  fit  to  the  values  of  2A^,  which  yields  a  slope  of  approximately  6.3. 

Examining  the  perturbation  to  the  other  eigenvalue,  -2A„  or  v=+0. 1374,  we  find 

only  one  difference.  The  summation  over  //  in  Eq  2.28  yields  so  there  is  an 

overall  difference  in  sign.  The  perturbed  eigenvalues  still  have  the  same  magnitude  as 
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the  perturbed  +2A^,  eigenvalue,  but  are  opposite  in  sign.  Thus,  linearizing  in  the 
vicinity  of  the  plus  branch  yields  eigenvalues  which  tend  toward  zero  with  increasing 
1//7. 

We  can  now  repeat  the  entire  process  which  led  to  Eq.  2.28-2.30,  but  perturbing 
about  the  minus  branch.  The  final  result  of  this  procedure  is  a  difference  in  sign  for  the 
perturbation.  Therefore,  linearizing  in  the  vicinity  of  the  minus  branch  yields 
eigenvalues  which  increase  in  absolute  value,  i.e.,  tend  away  from  zero  with  increasing 
1//T. 

Reaction  of  the  atom  upon  the  field  causes  several  other  effects  which  can  be 
explored  in  the  same  manner.  These  include  generation  of  extra  frequency 
components  in  the  field,  alteration  of  the  magnitudes  of  the  frequency  components, 
and  so  forth.  Several  of  these  have  been  explored  and  met  with  results  similar  to  those 
above.  Therefore,  they  need  not  be  discussed  further. 

We  are  now  prepared  to  return  to  the  question  of  chaos  in  the  semiclassical 
system.  As  the  Floquet  exponent  changes,  it  can  intersect  with  other  frequencies  in 
the  system.  For  the  plus  branch,  the  peaks  merge  near  «=90,  as  shown  in  Fig.  10.  As 
the  peaks  merge,  the  system  becomes  chaotic.  The  “wrong”  parity  frequency 
components  result  from  the  merging  of  the  different  frequency  components. 
Components  which  previously  had  no  means  of  interaction  can  now  interact  and  the 
result,  literally,  is  chaos. 

Figure  1 1  compares  the  Euclidean  separation  of  nearby  trajectories  for  7=85 

and  7=95.  For  7=95,  the  trajectory  is  clearly  non-chaotic — the  trajectories  actually 

become  closer  together  with  time,  then  return  to  their  original  separation.  On  the 

semi-logarithmic  scale,  there  are  no  linear  regions  indicating  exponential  growth. 

Conversely,  the  7=85  trajectory  exhibits  behavior  indicating  it  is  on  the  verge  of 

becoming  chaotic.  The  separation  exhibits  exponential  growth  for  approximately  1.5 

decades  before  decreasing.  Although  the  growth  does  not  remain  exponential  for 
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long,  the  initial  separation  is  definitely  exponential  and  thus  satisfies  the  definition  of 
chaos.  In  fact,  a  fit  to  the  linear  portion  of  the  graph  gives  an  exponent  0.042. 

However,  we  can  consider  the  other  indicators  of  chaos  as  well.  Figure  12 
shows  the  sequence  of  estimated  Lyapunov  exponents.  Like  the  Euclidean  distance, 
we  see  regions  of  exponential  growth  punctuated  by  exponential  decay.  In  order  to 
attempt  to  understand  this  phenomenon,  the  time  series  of  x  is  plotted  on  the  same 
graph.  The  time  series  shows  the  regions  of  exponential  growth  occur  when  the 
trajectory  is  close  to  the  Autler-Townes  solution.  As  the  trajectory  deviates,  as 
evidenced  by  the  increase  in  the  average  value  of  the  oscillations,  the  growth  is 
damped,  only  to  return  as  the  trajectory  returns  to  the  vicinity  of  the  Autler-Townes 
solution.  For  this  reason  we  claim  the  frequency  shift  from  the  Autler-Townes 
solution  is  related  to  chaos  in  the  following  way:  The  Autler-Townes  trajectory  is  a 
periodic  trajectory  which  becomes  unstable  as  a  result  of  these  shifts  in  the  frequency. 
Another  indication  of  chaos  is  the  qualitative  change  in  the  spectrum  given  by  the  FFT 
of  the  time  series,  as  shown  in  Fig.  13.  Observe  the  contrast  between  the  spectrum  for 
«=95  and  «=85.  Like  the  spectrum  shown  earlier  for  «  =  1024,  the  spectrum  for 
«  =95  consists  of  a  few  well  separated  peaks  which  can  be  identified  with  various 
harmonics  of  co'  and  2A„.  The  spectrum  for  /7=85,  on  the  other  hand,  shows  many 
more  peaks,  including  peaks  of  the  “wrong”  parity,  as  discussed  earlier.  Like  the 
Euclidean  separation,  the  FFT  does  not  show  the  system  to  be  strongly  chaotic. 
However,  it  does  show  a  qualitative  change  has  occurred  and  the  transition  to  chaotic 
behavior  is  underway. 

Finally,  consider  the  Poincare  section.  For  It  =95,  the  Poincare  section  is  a  very 
small  oval,  as  seen  in  Fig.  14a.  The  scatter  in  the  points  is  due  to  the  limited  precision 
in  interpolating  to  find  the  locations  of  the  zeroes  of  x(i).  Notice  how  small  an  area  in 
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x(t) 

Lyapunov  exponents  x  1 00 


Figure  12.  Sequence  of  Lyapunov  exponents  for  £=1,  «=85,  starting  near  the 
plus  branch,  (a)  Dotted  lines  bracket  possible  values  for  the  Lyapunov  exponent;  the 
dashed  line  indicates  a  linear  fit.  (b)  Overlay  of  x(i)  with  Lyapunov  exponent,  showing 
exponential  growth  in  separation  occurs  when  the  trajectory  is  near  Autler-Townes 
_ solution  (i.e.  where  oscillations  are  symmetric  about  r=0). 
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the  phase  space  is  enclosed.  However,  Fig,  14b-d  shows  that  as  Ji  decreases  the  area 
enclosed  by  the  Poincare  section  increases  in  size  and  complexity. 

All  three  tests  for  chaos  show  a  qualitative  change  as  77  decreases  from  95  to  85. 
This  transition  is  marked  by  the  merging  of  two  frequency  components  which  results 
in  the  system  becoming  chaotic.  The  frequency  shifts  can  be  traced  to  an  interaction 
predicted  by  second-order  perturbation  theory.  Of  course,  a  single  instance  would  not 
be  sufficient  to  establish  a  cause  and  effect  relationship,  so  let  us  now  consider  the 
transition  to  chaotic  behavior  for  the  minus  branch  with  the  same  value  of  s.  The 
magnitude  of  the  slope  of  the  minus  branch  line  in  Fig.  9b  is  approximately  6.5,  in 
reasonable  agreement  with  the  value  given  by  Eq.  2.30.  There  are  no  frequency 
components  appearing  in  Eq,  2.9  which  are  slightly  larger  than  2co'-2X^.  However,  we 
do  observe  in  the  spectrum  of  x  a  harmonic  at  -(3<y -4/?.^)=0.7272.  These  two  peaks 
merge  near  77=32,  and  the  transition  to  chaos  is  much  stronger.  Using  the  same  tests 
for  chaos,  let  us  examine  this  transition.  First,  the  Euclidean  separation  and  Lyapunov 
exponent  are  shown  in  Fig.  1 5.  Just  as  for  the  plus  branch,  we  observe  a  change  in  the 
behavior  of  the  system  as  77  drops  from  35  to  30.  The  Euclidean  separation  shows  an 
exponential  growth,  which  in  this  case  does  not  entirely  decay  before  the  next  growth 
cycle  occurs.  An  exponential  fit  yields  an  exponent  of  approximately  0.050.  The 
separation  between  non-chaotic  trajectories,  while  growing,  is  not  exponential.  This  is 
very  similar  to  the  behavior  seen  in  the  Henon-Heiles  system  in  Sect.  2. 1 .  The 
Lyapunov  exponent  for  the  chaotic  trajectory  converges  to  approximately  0.015.  The 
rapidity  with  which  the  Lyapunov  exponent  converges  shows  the  more  robust  nature 
of  the  chaos  in  this  portion  (smaller  77)  of  the  parameter  space. 


46 


47 


48 


The  difference  in  the  FFTs  for  77=30  and  77=32  also  demonstrates  the  change 
from  regular  to  chaotic  behavior.  The  two  peaks  at  <y=0.304  and  <^=0.378  are  the  two 
harmonics  which  merge  (3a>'-4A  and  2(v'-2A^)  as  discussed  above.  They  have  shifted 
significantly  from  their  77=  c»  values  of  <y=0.7272  and  £y=0. 1364.  As  with  the  plus 
branch,  the  shift  departs  from  linearity  shortly  before  the  peaks  merge.  As  they  do,  the 
frequency  components  become  broadband  regions  instead  of  individual  peaks, 
demonstrated  by  Fig.  16.  Incidentally,  the  “wrong”  parity  frequency  components  do 
not  become  significant,  probably  because  the  harmonics  did  not  merge  in  the  vicinity 
of  the  frequency  components  which  were  absent.  Finally,  we  observe  the  change  in 
the  Poincare  section  in  Fig.  17.  In  this  instance,  the  break  from  regular  to  chaotic 
behavior  can  be  observed  with  increasing  time.  For  77=30.2,  the  transition  occurs 
near  /  =  621 .  Before  this  time,  all  the  points  in  the  Poincare  section  fall  along  the 
curve  formed  by  the  +  marks.  Afterwards,  the  points  are  more  scattered,  falling 
generally  in  the  same  region  as  for  the  plus  branch. 

The  reader  should  notice  an  important  difference  exists  between  the  transition  to 
chaos  in  the  plus  and  minus  branches.  Namely,  in  the  minus  branch,  one  of  the 
merging  harmonics  does  not  exist  in  the  Autler-Townes  model.  This  harmonic,  like 
the  other  higher  harmonics  of2A^  in  Table  1,  result  from  the  term  which  was  discarded 
in  Eq.  2. 17.  The  relatively  low  value  of  77  suggests  that  the  Autler-Townes  model  is 
of  relatively  little  relevance;  the  success  of  the  approach  here  is  not  so  much  the 
prediction  of  chaos,  but  rather  the  prediction  that  chaos  occurs  sooner  for  the  plus 
branch  than  for  the  minus  branch. 

Our  analytical  approach  has  three  important  limitations.  First,  the  Autler- 

Townes  trajectory  is  not  identical  to  the  periodic  trajectory  in  the  semiclassical  model. 

It  is  close  to  the  periodic  trajectory,  but  the  starting  points  are  not  exactly  the  same. 

especially  for  smaller  77.  Second,  perturbation  theory  gives  only  the  lowest  order 

correction  in  1/77,  i.e.,  the  slope  of  the  1/77  plot.  Third,  we  have  neglected  an 
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inhomogeneous  term  which  generates  harmonics  which  this  approach  is  unable  to 
account  for.  Nonetheless,  these  results  do  demonstrate  the  transition  to  chaos  occurs 
by  reducing  n  from  77=oo,  where  the  semiclassical  model  is  identical  to  the  Autler- 
Townes  model,  to  i7<100,  where  chaos  begins  to  emerge.  Furthermore,  we  are  able 
to  attain  a  qualitative  understanding  of  why  one  branch  becomes  unstable  before  the 
other.  Before  moving  ahead  to  examine  the  relationship  between  chaos  as  it  has  been 
observed  here  and  the  effects  on  the  full  quantum  model,  we  must  examine  the 
relationship  between  the  semiclassical  model  and  the  quantum  model  in  the  non- 
chaotic  region  of  parameter  space. 
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CHAPTER  3 


QUASICLASSICAL  STATES  AND  THE  QUANTUM 

MODEL 

In  Chapter  2  we  pointed  out  that  for  a  semiclassical  treatment  to  be  valid,  the 
atomic  and  field  portions  of  the  wave  function  must  factorize.  Assuming  that  they 
factorize  and  then  treating  the  expectation  values  of  the  field  variables  as  dynamical 
variables  makes  the  approximation  semiclassical.  States  which  actually  do  remain 
approximately  factorizable  have  been  termed  quasiclassical.  A  quasiclassical  state  will 
have  a  relatively  well  defined  field  amplitude  and  phase.  If  A/?  is  the  uncertainty  in  the 
photon  number  of  the  state  and  is  the  uncertainty  in  the  field  phase,  then  we  require 

A///«  «1  and  A^<1 .  These  conditions  can  be  satisfied  by  coherent  states  with 
sufficiently  large  n  and  squeezed  states  in  which  the  squeezing  is  not  too  large.  In  this 
work  we  will  consider  only  initial  conditions  which  involve  coherent  states.  However, 
the  fact  that  the  approximation  remains  good  for  a  range  of  times,  even  though  the 
field  does  not  remain  in  a  coherent  state  after  /=0,  demonstrates  the  broader 
applicability  of  the  technique. 

In  Sec.  1  of  this  chapter  we  will  first  give  a  general  development  of 
quasiclassical  states  and  show  how  the  plus  and  minus  branches  of  the  Autler-Townes 
solutions  relate  to  the  quasiclassical  states.  In  Section  2  we  will  examine  the  evolution 
of  the  field  and  show  how  the  plus  and  minus  branches  may  be  used  to  obtain  the 
evolution  for  arbitrary  initial  states.  Section  3  will  give  an  alternate  development  using 
just  the  eigenstates  of  the  quantum  model. 
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3.1  Development  of  Quasiclassical  States 


Let  a  quantum  harmonic  oscillator  interact  with  an  atomic  system.  If  the 
oscillator  is  sufficiently  energetic,  the  atom  cannot  strongly  influence  the  total  energy 
of  the  oscillator.  As  mentioned  in  Chapter  2,  we  desire  to  find  states  such  that 
I  ^(0))  =  I  y/(0))  0 1  a)  ->  I  T(/))  =  I  y/(/))  ®  1 0(/)),  (3.1) 

(where  1T(/))  is  the  total  wavefunction,  |  y4j))  is  the  atomic  part  and  |0(/))  is  the  field 

part)  for  times  />0.  The  atomic  system,  which  could  be  a  two-  or  three-level  system 
or  even  a  few  such  systems,  has  a  state  space  of  dimension  N.  There  can  then  be  no 
more  than  Af  quasiclassical  trajectories,  as  follows;  Given  Af  quasiclassical  trajectories, 
any  other  initial  condition  can  be  written  as  a  linear  combination  of  the  // linearly 
independent  (though  perhaps  not  orthogonal)  trajectories.  The  new  initial  condition 
will  evolve  according  to  the  linear  combination  of  N  wavefunctions  of  the  form  (3.1). 
This  linear  combination  will  not  be  generally  factorizable. 

Classically,  if  two  interacting  systems,  such  as  the  atom  and  the  field,  exchange 
little  energy,  we  expect  the  smaller  system  to  have  little  effect  on  the  larger  system. 

But  quantum  mechanically  the  systems  can  interact  in  a  manner  which  need  not  involve 
substantial  energy  exchange,  through  entanglement.  By  becoming  entangled,  the 
smaller  system  can  exert  significant  influence  on  the  larger  system  by,  for  example, 
causing  the  larger  system  to  evolve  into  Schrodinger  cat  states.  By  avoiding 
entanglement  (i.e.  remaining  pure  or  factorizable),  quasiclassical  states  avoid  this 
quantum  mechanical  phenomenon  and  thus  the  oscillator  will  be  least  affected.  We 
expect  the  oscillations  to  be  approximately  constant,  but  we  will  allow  for  a  slight 
modification  of  the  frequency  co  to  <y'  (as  we  observed  in  Chap.  2).  In  order  for  the 
system  to  remain  factorizable,  we  will  require  that  at  the  end  of  each  period  of  the 
oscillator,  the  atomic  system  must  return  to  very  near  its  initial  state,  up  to  an  overall 
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phase  factor.  To  put  it  another  way, 

I  v/.(/  =  2nlo}'))  =  y/,(0)>.  (3.2) 

Taking  the  form  of  Eq.  3.2  along  with  the  assumption  the  oscillator  is  not 
significantly  disturbed  leads  naturally  to  the  application  of  Floquet’s  theorem  to  find 
the  quasiclassical  states.  Making  the  approximation  that  the  oscillator  is  unaffected 
(other  than  the  shift  in  frequency)  leads  to  the  prescribed  field  Hamiltonian 

H  =  +2fig^Jno-^  cosco'/ ,  (3.3) 

which  is  the  Autler-Townes  Hamiltonian  introduced  in  Sect  1.1.  Now  the  connection 
between  the  Autler-Townes  states  and  the  quasiclassical  states  is  clear:  the 
quasiclassical  states  must  have  Floquet-like  behavior,  while  the  Autler-Townes  states 
exactly  solve  a  Floquet  Hamiltonian.  A  comparison  between  the  Autler-Townes 
trajectories  for  «=32,£=0. 141  (g=0.05)  and  the  numerical  integration  of  the  fiilly 
quantized  system  is  shown  in  Fig.  18.  The  atom  starts  on  the  Autler-Townes  plus 
branch,  while  in  this  figure  and  elsewhere  unless  noted  otherwise,  the  fully  quantized 
field  starts  in  a  coherent  state  with  a  =  .  The  interaction  picture,  in  which  the  axis 

system  of  the  Bloch  sphere  rotates  at  the  frequency  (Vq,  is  useful  for  small  fto  observe 
the  interaction-induced  detuning.  The  azimuth  angle  traversed  by  a  particular  point 
along  the  teardrop-shaped  trajectory  gives  the  total  drift  {(o'-co)t.  For  large  e,  the 
trajectories  become  complicated  (recall  Fig.  6)  and  the  interaction  picture  is  not 
particularly  useful.  In  this  case,  though,  the  interaction  picture  demonstrates  the 
similarity  of  the  detuned  Autler-Townes  model  and  the  fully  quantized  model, 
particularly  once  the  quantum  model  settles  down.  The  total  time  shown  is 
approximately  P=2Q0  (in  units  with  co  scaled  to  1)  or  about  32  cycles.  The  first  5 
cycles  or  so  show  some  differences  between  the  models,  but  after  that  time  they  are 
virtually  identical.  The  agreement  is  remarkable,  considering  the  difference  in 
computations  required:  the  quantum  model  requires  integration  of  nearly 
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(a) 


(b) 


Figure  1 8.  Comparison  of  Autler-Townes  model  with  detuning  and  quantum  model. 
In  all  plots  £=0.141,  «=32,  g=0.05.  (a)  Autler-Townes  model,  Schrodinger  picture, 
(b)  Autler-Townes  model,  interaction  picture  (rotating  with  angular  frequency  oS).  (c) 
_ and  (d)  are  the  results  for  the  fully  quantized  model. 


250  equations,  while  the  Autler-Townes  results  can  be  obtained  by  numerically 
integrating  3  equations  or  using  the  closed  form  expression  given  by  Eq.  2.9,  together 
with  the  y4’s  and  B’s  from  Ref  9.  Agreement  improves  with  larger  77,  where  the 
quantum  model  requires  even  more  effort  (roughly  4-(77  +  5-\/^)  equations). 

The  factorization  into  an  atomic  part  and  a  field  part  holds  for  a  fairly  long  time, 
as  shown  by  the  atomic  state  purity  in  Fig.  19.  The  purity  is  measured  by  the  T^(p^^^), 
where  Pat  =  Trr,e|(j(/?)  and  p  is  the  density  operator.  A  purity  of  1  indicates  a  pure 
state,  while  the  minimum  purity  for  a  two  level  system  is  0.5  and  indicates  a  maximally 
entangled  system.  In  Fig.  19  the  purity  remains  above  0.98  until  about  7=400.  This 
time  is  about  twice  the  time  shown  in  Fig.  1 8.  To  fully  appreciate  the  value  of  the 
quasiclassical  states,  consider  the  trajectories  generated  by  starting  in  the  state  \a) 

instead  of  1^+).  Figure  20a,b  shows  the  Schrodinger  picture  evolution  of  the 

semiclassical  system  (2.5)  and  the  fully  quantized  system  for  the  same  parameters  as  in 
Fig.  18,  but  with  the  atom  initially  in  the  upper  state.  The  difference  in  the  two 
trajectories  is  immediately  obvious  as  the  well-known  “collapse”  of  the  wavefunction 
occurs.2^>38.40^2  Yhe  collapse  is  even  more  apparent  in  Fig.  20c,  where  the  purity 
rapidly  drops  off  The  approximation  of  factorization,  therefore,  is  seen  to 
dramatically  fail  for  states  which  are  not  near  the  quasiclassical  states.  The  collapse  of 
a  general  state  is  a  result  of  the  field  splitting  into  two  parts,'*^’'^  one  associated  with 
the  plus  branch  and  one  associated  with  the  minus  branch.  Figure  2 1  shows  the 
quasiprobability  Q  distribution  at  the  beginning  and  end  of  the  trajectories  shown  in 
Fig.  20.  The  Q  distribution  is  defined  as 

0(«)  =  («ke.d|«)  (3.4) 

where  is  obtained  from  the  total  density  operator,  p,  by  tracing  over  the  atomic 
states.  The  variable  a  represents  the  complex  amplitude  of  a  coherent  state.  This 
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figure  is  in  the  interaction  picture,  rotating  with  angular  frequency  -co.  Thus  the  plus 
branch  portion  of  the  field,  which  rotates  with  co'-o}>0  is  the  lower  peak  (rotating 
clockwise),  while  the  minus  branch  portion  of  the  field  is  the  upper  peak  (rotating 
counter-clockwise).  Referring  to  Fig.  7,  one  can  see  the  state  \a),  which  has  is 

closer  to  the  minus  branch  than  the  plus  branch.  As  a  result,  the  minus  branch  portion 
of  the  field  has  a  greater  weight,  and  so  the  peak  is  bigger  (has  more  contours)  in  Fig. 
21 .  This  association  is  demonstrated  in  Fig.  22,  where  the  atom  starts  on  the  plus 
branch  and  only  one  peak  is  observed,  rotating  clockwise  in  the  interaction  picture. 

We  will  discuss  this  further  in  Sec.  3.2. 

The  location  of  the  peak  of  the  Q  distribution  approximately  gives  the  value  of 
the  complex  field  phasor  discussed  at  the  end  of  Sec.  2.2.  Using  Fig.  22,  a^3.76- 
4.32/.  This  yields  a  detuning  of  0. 855/20 l.l»i  0.0043  (please  recall  in  numerical 
calculations,  time  has  been  scaled  so  that  oy=\).  The  value  for  the  interaction  induced 
detuning  given  by  Eq.  2. 14  is  0.00427,  in  complete  agreement  at  this  level  of 
precision.  Greater  precision  can  be  obtained  by  calculating  the  expectation  values  of 
the  operators  and  the  real  and  imaginary  parts  of  the  boson  annihilation 
operator,  and  using  the  same  procedure  for  constructing  the  complex  field  phasor. 

This  technique  yields  a  value  for  the  detuning  of  0.00422,  which  is  in  remarkable 
agreement  considering  the  differences  in  the  models.  In  fact,  this  technique  provides 
good  agreement  over  a  fairly  broad  range  of  values  for  e.  Letting  the  product  of 
thedetuning  and  the  average  photon  number  be  given  by  57,  Figure  23  plots  57  versus 
^The  solid  line  gives  the  value  from  Eq.  2. 14,  while  the  squares  and  circles  give  the 
values  calculated  from  the  quantum  model  for  77=32  and  64,  respectively.  This  figure 
demonstrates  the  improved  agreement  which  can  be  obtained  by  increasing  77 . 
Incidentally,  if  this  procedure  was  used  to  calculate  the  detuning  in  the  quantum  model 
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Figure  23.  Comparison  of  the  interaction  induced  detuning  from  the 
semiclassical  model  (Eq.  2.14)  and  the  fully  quantum  model  along  the  plus 
branch.  Agreement  improves  with  increase  n . 
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for  a  situation  like  that  in  Fig.  20,  the  value  would  be  almost  zero,  since  the 
contributions  from  the  two  peaks  would  nearly  cancel  each  other. 

Like  the  correspondence  between  the  Autler-Townes  model  and  the 
semiclassical  model,  the  correspondence  between  the  quasiclassical  states  and  the 
semiclassical  model  improves  with  increasing  n .  Using  Eq.  2.4,  the  expansion  ofx 
and  as  complex  exponentials,  and  Ref  38  for  the  series  expansion  of  the  Autler- 
Townes  A  and  B  coefficients  in  s,  we  find  that 


(a)=3Q)') 

a\ 


(3.5) 


Therefore,  the  third  and  higher  harmonics  can  be  made  negligible  given  sufficiently 
large  ti .  Thus,  the  agreement  between  the  quantum  solution  along  a  particular  branch 
and  the  Autler-Townes  Floquet  eigenstates  (or  the  semiclassical  solution)  can  be  made 
as  close  as  desired  for  fixed  shy  increasing  7i .  On  the  other  hand,  by  fixing  the 
agreement  deteriorates  with  increasing  s. 


3.2  Evolution  of  the  Field  and  of  Arbitrary  Initial  States 


Up  to  this  point,  almost  all  the  attention  has  been  focused  upon  the  evolution  of 
the  atom.  Now  we  shall  consider  the  evolution  of  the  field  for  the  quantum  model. 
The  semiclassical  model  gives  expectation  values  for  the  field  quadratures  and  ^2 
along  a  single  branch;  however,  we  need  more  information  about  the  quantum  field 
state  in  order  to  calculate  the  evolution  of  linear  combinations  of  the  two  branches. 

To  find  how  the  quantum  field  evolves,  assume  a  factorized  form 

|HK^(/))  =1  ^^',.(/))0|O^(/))  for  the  total  wavefunction  and  insert  it  into  the  time- 

dependent  Schrodinger  Equation: 

^(1  n  (0)  (0))  =  -^  //(I  K'.  (0)  ®  I  <!>.  (' ))) .  (3.6) 
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Now  project  out  the  portion  which  follows  along  the  atomic  quasiclassical  trajectory 

=  -j(»'.(')|//k.(0)|®,(0> 

On  the  RHS  of  the  first  line  of  Eq.  3.7,  the  derivative  w.r.t.  time  of  the  atomic  state 
appears.  This  derivative  is  given  by  the  action  of  the  Floquet  Hamiltonian  (including 
the  interaction  induced  detuning) 

f(K'.(')|||v'.(0)j  =  -^(r.(')i//.,|^.(')>.  (3.8) 

Making  this  substitution  in  Eq.  3.7,  we  obtain  an  equation  for  the  evolution  of  the  field 

||<l>.(0)  =  -|(K'.(0|(«  -  1<'.(0>|1>.(')>  (3.9) 

and  then  substitute  H  and  from  Eq.  1.1  and  Eq.  1 .2  to  get 

= -/  {txi^a  + gx^[i){a-\-a^ cosa)'i)\Q>^{t)).  (3.10) 

In  this  manner  the  technique  has  come  full  circle:  We  started  by  assuming 
factorization  and  found  the  atomic  evolution  using  the  expectation  values  for  the  field 
in  a  coherent  state  \a)  with  a  =  Then  we  used  the  expectation  values  for  the 

atom  and  found  the  evolution  of  the  field. 

The  full  solution  to  Eq.  3. 10  is  a  coherent  state  times  a  phase  factor  and  is  given 
in  Ref  38.  Here  we  note  only  that  to  lowest  order  in  s/ji,  the  complex  amplitude  of 
the  coherent  state  is 

=  l  +  O  .  (3.11) 

Similar  results  are  found  for  the  minus  branch  by  replacing  the  plus  branch  expectation 
values  with  those  for  the  minus  branch  and  using  the  value  of  co'  for  the  minus  branch. 
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Similar  analysis  for  quasiclassical  states  using  the  rotating  wave  approximation 
(RWA)  has  shown  the  field  states  do  not  remain  coherent  for  very  long  times,  but 
rather  exhibit  squeezing.  jn  fact,  inspection  of  Fig.  22  shows  the  state  is  no  longer 
coherent,  but  has  increased  phase  uncertainty.  However,  numerical  calculations  show 
the  coherent  state  approximation  to  be  adequate  for  short  enough  times  (see  below). 

As  discussed  at  the  beginning  of  this  chapter,  the  existence  of  quasiclassical  states  does 
not  depend  on  the  precise  form  of  the  field,  so  long  as  the  phase  and  amplitude  remain 
fairly  well  defined. 

Having  an  expression  for  the  field  evolution  along  each  branch,  we  can  now 
move  on  to  the  evolution  of  the  field  and  atom  for  an  arbitrary  initial  atomic  state. 

Since  the  two  branches  at  1=0  are  orthogonal,  any  atomic  state  can  be  written  as  a 
linear  superposition  of  the  two  branches.  Using  the  same  form  as  in  Eq.  2.9,  we  can 
write 


I'F(o))  =  [ cos-^l  y/S0))  +  e‘* s\n^\  y/_(o) 


(3.12) 


The  detuning  for  each  branch  is  opposite  in  sign,  so  the  field  states  drift  apart  until 
they  become  approximately  orthogonal.  In  terms  of  the  Q  distribution,  the  states  are 
orthogonal  when  the  two  peaks  are  well  separated,  i.e.,  they  do  not  overlap.  The 
width  of  the  distribution  Q{a)  for  a  coherent  state  |cif')  is  independent  of  the  complex 


amplitude  a\  and  is  approximately  equal  to  unity.  However,  from  the  origin,  the  angle 
subtended  by  the  distribution  (i.e.,  the  phase  uncertainty)  falls  as  l/|cr'|  =  \lyfn  . 

Therefore,  the  two  branches  will  be  roughly  separated  at  a  time  such  that 
|ry'-ry|/^  =  UseofEq.  2.14  for  the  detuning  yields 


oM 


c 


(3.13) 
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where  Xq  is  amplitude  of  the  first  term  in  the  cosine  expansion  (2. 12)  of  x{/)  along  one 
of  the  branches.  This  time  is  the  collapse  time,  referring  to  the  collapse  of  the 
waveflinction  mentioned  at  the  beginning  of  Sect  3.1.  It  is  called  a  collapse  because  in 
the  RWA,  after  the  collapse  the  Rabi  oscillations  in  the  inversion  cease.  However,  the 
collapse  is  actually  the  collapse  of  the  interference  terms  between  the  field  states  along 
each  branch.  For  a  fixed  value  of  £,  which  fixes  the  value  of  Xq,  the  time  to  achieve 
collapse  grows  with  77.  For  the  case  shown  in  Fig.  18,  Eq.  3.13  yields  <y/^=41.2, 
which  agrees  well  with  the  time  of  the  collapse  shown  in  Fig.  20.  Alternatively,  for 
fixed  g,  the  collapse  time  depends  only  indirectly  upon  77  (through  Xq).  This  result  is 
comparable  to  the  result  for  the  RWA,  for  which  is  independent  of  77.  After  the 
collapse  occurs,  the  expectation  values  of  most  quantities  are  given  by  a  weighted 
average  along  each  branch  separately  (with  weights  cos^y/l  and  siv?-}^2  for  the  plus 
and  minus  branches,  respectively). 

Finally,  the  two  field  states  will  overlap  as  \(o'  -  <y|/^  =  ti.  This  overlap  will 
result  in  a  revival  of  interference,  which  will  occur  at 


col„  = 


mi 


iKyfn 


(3.14) 


This  result  agrees  with  the  work  done  by  Zaheer  and  Zubaiiy**^.  For  g=\/ sVTo , 
77=10,  £=0.1,  they  show  in  their  Fig.  2  a  revival  which  peaks  near  2gt=40,  or  cot=3l6. 
Equation  3.14  yields  a  value  of  <y/y^=3 1 9.  However,  this  small  value  of  e  is  insufficient 
to  distinguish  these  results  from  the  RWA,  for  which  the  revival  time  is  given  by  Eq. 
3.14  with  Xq  set  equal  to  unity.  The  RWA  yields  for  these  parameters  <y7y^=3 14.  Using 
a  larger  value  of  £=0.3  and  77=36,  Fig.  24  shows  a  rather  indistinct  revival  occurring  in 
the  vicinity  of  <y/=500.  For  these  values,  Eq.  3.14  yields  a  revival  time  of  cotff=49\, 
while  the  RWA  yields  (otj^311,  in  sharp  contrast.  In  order  to  better  determine  when 
the  revival  occurs.  Fig.  25  shows  the  contribution  to  the  inversion  from  just  the 
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interference  terms.  The  interference  contributions  can  be  found  by  performing  three 
separate  numerical  integrations  of  the  quantum  model.  Let  (z+)  be  the  expectation 
value  given  by  the  evolution  of  )y/+(/)),  (z,)  from  |y/(/)),  and  (z)  from  the  linear 


combination  of  the  two  (as  already  shown  in  Fig.  24).  The  interference  terms  will  then 
be  equal  to  : 


(3.15) 


This  quantity  is  plotted  in  Fig.  25,  and  the  revival  of  the  interference  between  the  two 
branches  is  clearly  evident,  becoming  largest  very  close  to  /=500.  Now  it  is  apparent 
that  the  results  coming  from  the  Floquet  eigenstates  agree  well  with  the  full  quantum 
model. 


The  indefinite  revival  demonstrates  that  the  coherent  state  approximation  is 
beginning  to  fail  for  this  large  time.  If  the  field  remained  in  a  coherent  state,  the 
revival  would  last  as  short  a  time  as  the  collapse,  but  Fig.  25  shows  that  it  happens 
over  a  substantially  longer  time.  Figure  26  shows  how  the  0  distributions  for  the  plus 
and  minus  branches  individually  evolve.  There  is  little  change  in  the  radial  profile, 
indicating  the  amplitude  of  the  coherent  state  and  its  uncertainty  do  not  change  much. 
However,  the  distribution  becomes  very  elongated,  corresponding  to  an  increased 
uncertainty  in  the  phase.  This  is  suggestive  of  the  squeezing  seen  in  the  RWA.'^^ 
Finally,  Fig.  27  shows  the  evolution  of  the  Q  distribution  for  a  linear  combination  of 
plus  and  minus  branches.  Comparison  of  Figs.  26  and  27  shows  the  Q  distribution  for 
the  linear  combination  evolves  very  much  like  the  sum  of  the  Q  distributions  for  the 
two  branches,  which  further  reinforces  the  picture  of  the  quasiclassical  states. 

An  attempt  to  improve  the  accuracy  of  the  coherent  state  approximation  at  times 
near  the  revival  time  by  increasing  7/  cannot  be  successful:  the  revival  time  increases 
linearly  with  77  for  fixed  s,  while  the  discrepancy  with  the  quasiclassical  states  (given 
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by  Eq.  3.5)  decreases  only  as  l/V^.  Although  the  discrepancy  given  by  Eq.  3.5  is  not 
necessarily  directly  connected  to  the  failure  of  the  coherent  state  approximation,  it 
gives  an  upper  bound  to  the  agreement  which  can  be  expected  between  the 
quasiclassical  states  and  the  full  quantum  model. 

Having  considered  some  of  the  limitations  of  the  use  of  the  quasiclassical  states, 
now  let  us  explore  a  concrete  example  of  their  use.  For  simplicity,  consider  the  small  s 
limit,  so  that  the  series  expansion  of/l’s,  5’s,  and  probability  amplitudes  (all  given  in 
Ref  38)  can  be  truncated  to  first  order  in  e.  A  particularly  useful  state  for 
consideration  is  the  upper  atomic  energy  eigenstate,  |a).  It  has  been  widely  used  as  an 

initial  condition  by  various  authors,  and  the  factorization  approximation  used  to 
develop  the  semiclassical  model  fails  dramatically  from  this  state,  as  already  seen. 
Following  Ref  38,  the  state  \a)  is  formed  by  the  linear  combination 


(3.16) 


Now  use  the  probability  amplitudes  in  Eq.  3.16,  along  with  Eq.  1 .3  for  the  plus  and 
minus  branches,  to  form 

=  i{(l  - +  [(l  -  (0)  (317) 

+ee'‘“-')\a)^-[\  +  e-  (!).(/)) 

with  co^  the  detuned  frequency  for  the  appropriate  branch,  which  is 
g  e 

=0}±  r-^  =  o}±  — 

2V^  li 

to  first  order  in  e.  Although  the  form  of  |d)+(/))  is  somewhat  complicated,  for  times 
smaller  than  the  revival  time  and  £//7  sufficiently  small,  it  simplifies  to  the  following 
product  of  coherent  state  and  phase  factor 

(3.18) 
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Typically  plots  of  the  inversion,  z(/).  are  shown  in  the  literature.  To  obtain  r(/), 
take  the  expectation  value  of  for  the  state  |'F(/)).  Making  these  substitutions  and 

performing  straightforward  algebraic  manipulation  yields  for  the  inversion  when 
starting  from  the  state  \a) 

z(/)  s  |^cos(2g-\/^/)  - 2Esm{2(ot) sin(2g^-v/^/)je"*’'’'' 

(  (3.19a) 

+2£‘sin(2<'<j^)sin 


or,  using  the  scaled  time  with  oy=\. 


z{t)  =[cos(4ff)-2£'sin(2/)sin(4f/)]e 


(3.19b) 


+2£’sin(2/)sin 

V  n 


The  first  form  shows  that  as  s-^0,  the  RWA  result  emerges;  Rabi  oscillations  at 
frequency  2g4^  modulated  by  a  Gaussian  collapse  envelope.  The  second  form  is 
useful  to  see  how  the  various  quantities  scale  with  fixed  fas  «  is  varied.  In  either 
equation,  the  interference  between  the  two  branches  is  responsible  for  the  terms  in 
square  brackets  (which  vanish  after  the  collapse  time).  This  result  does  not  show  any 
revivals  because  of  the  approximations  made  in  obtaining  Eq.  3.18.  However,  it  does 
show  the  zero  in  the  oscillations  at  one-half  the  revival  time,  which  can  be  observed  in 
Fig.  24.  A  comparison  between  the  results  of  numerically  integrating  the  full  quantum 
model  and  Eq.  3.19b  is  shown  in  Fig.  28.  The  trajectories  in  this  figure  use  the  same 
parameters  as  in  Fig.  18;  £=0. 141,  ^0.05,  n=32.  Even  for  a  value  of  f  this  large,  the 
agreement  is  quite  good. 

So  far  in  this  chapter  we  have  seen  how  the  Floquet  eigenstates  which  come 
from  the  Autler-Townes  model  are  connected  to  quasiclassical  states  which  remain 
factorizable  for  long  times.  We  then  found  that  any  state  could  be  written  as  a  linear 
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combination  of  the  plus  and  minus  branches  and  its  evolution  could  then  be  calculated 
in  a  straightforward  manner  for  times  which  approach  the  revival  time.  However, 
there  is  an  additional  manner  in  which  these  special  states  can  be  derived  for  the 
quantum  system,  without  any  reference  to  the  semiclassical  system. 

3.3  Quasiclassical  States  From  Eigenvectors  of  the  Hamiltonian 


The  motivation  for  finding  the  quasiclassical  states  developed  from  a  need  to 
find  states  which  closely  matched  the  evolution  of  the  semiclassical  system.  However, 
it  is  also  possible  to  find  such  states  using  an  approximate  treatment  of  the 
eigenvectors  of  the  Hamiltonian  in  Eq.  1.1.  Graham  and  Hohnerbach  published  a 
solution  for  the  eigenvectors  to  zeroth  order  and  the  eigenvalues  correct  to  first 
order.26  Here  we  present  a  similar,  though  somewhat  different  technique  to  obtain  the 
eigenvalues  and  eigenvectors  correct  to  first  order. 

The  term  to  be  treated  perturbatively  is  the  atomic  energy.  For  field  states 
which  are  built  primarily  from  a  limited  range  of  field  eigenstates  |/7),  such  as  coherent 

states,  a  w  «  yfn .  If  the  coupling  and  the  average  photon  number  are  not  too  small, 
such  that  gyfn  >  j  cOq  or  £>j,  the  interaction  term  will  be  larger  than  the  atomic 

energy.  Dropping  the  atomic  energy  term  from  Eq.  1 . 1  we  obtain  the  unperturbed 
Hamiltonian 

=  ticoa^a  +  Tiga^{a+a^).  (3.20) 


Defining  the  atomic  states 


|±x)  =  ^(|<,)±|i)) 


(3.21) 


which  are  the  eigenstates  of  cr^  with  eigenvalue  ±1  (and  hence  lie  at  x=±I  on  the  Bloch 
sphere),  we  can  write  for  the  eigenstates  of  Hq 
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]ji(ixi^a^tigaX<3  +  <^^)  |±^)|^«±)=  ti(ai^a±Tig[a  +  a^)  |+x)|(I>„j) 


(3.22) 


with  the  field  state,  which  is  still  to  be  determined.  Now  replace  a  and  a+  with 

new  canonical  annihilation  and  creation  operators  a'=a±glo)  and  a'^=a^±gl(o,  with  the 
sign  chosen  to  match  the  sign  of  the  atomic  state.  Upon  substitution  the  eigenvalue 
equation  becomes,  (after  projecting  out  the  atomic  part) 


O  )-£  O  ) 

^ntf  ^n±/* 


(3.23) 


Another  way  to  write  the  new  annihilation  operator  is  a±'=D'^(±g/(jo)aD(±g/co),  where 
D(a)  is  the  displacement  operator  which  generates  the  coherent  state  \a)  when  applied 

to  the  vacuum  state.  Using  this  definition,  Eq.  3.23  can  be  written 


na}D^\±^ala,D\±^  -  Tio}\ 


\(0 


(3.24) 


Now  the  eigenstates  can  be  determined  by  inspection.  If  |/?)  is  a  number  state,  then 


(3.25) 


along  with  the  corresponding  eigenvalues 


=  ntico-tico 


(3.26) 


which  are  doubly  degenerate,  since  they  are  independent  of  the  sign  of  the  state.  First 
order  degenerate  perturbation  theory  requires  that  we  calculate  the  matrix  elements  of 
the  perturbation  within  each  degenerate  subspace.  Diagonalization  of  the  resulting 
matrix  yields  the  zeroth  order  eigenvectors  and  the  first  order  correction  to  the 
eigenvalues.47  In  this  instance, 
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(3.27) 


(+x|a,|+x) 

Using  Eq.  3.21,  it  is  straightforward  to  evaluate  the  atomic  state  matrix  elements.  The 
diagonal  elements  are  both  zero,  and  the  atomic  part  of  the  off-diagonal  elements  is 
equal  to  1.  The  field  part  is  more  difficult,  but  de  Oliveira,  e(  ai,  have  reported 
thoroughly  on  the  displaced  number  states."*^  The  form  of  the  displacement  operator 
Z3(a)=exp(<awr+-«*a)  immediately  shows  that  D^{a)=D{-a).  From  Ref  48,  then,  we 
obtain 


(3.28) 


where  are  the  Laguerre  polynomials.  Since  this  quantity  is  real  valued,  the  matrix 
elements  are  equal.  The  eigenvalues  correct  to  first  order  are  then 


E  ^  =  Tio) 

n±  I 


"-{if 


(3.29) 


and  the  zeroth  order  eigenvectors  are 


v2  L  v<y/  \  0) j 


(3.30) 


One  may  wonder  what  a  displaced  number  state  looks  like.  Reference  48  focuses  on 
states  with  a  large  displacement  and  small  photon  number,  while  the  opposite  is 
needed  here.  Briefly,  for  large  n  and  fairly  small  glo),  the  photon  number  probability 
distribution  is  narrow  and  centered  near  //.  As  glco  increases,  oscillations  in  the 
distribution  develop  and  it  becomes  wider,  as  demonstrated  in  Fig.  29.  Analytically, 
the  average  photon  number  is 


(/;)  =  /7  + 


(3.31) 
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Figure  30.  Eigenvalues  from  numerical  diagonalization  of  Hamiltonian  (“exact”),  Eq 
3.26  (zeroth  order),  and  Eq.  3.29  (first  order)  for  g/fo  of  (a)  0.05  and  (b)  0.25.  The 


values  plotted  are  vs.  //. 


and  the  uncertainty  in  the  average  photon  number  is 


((A»)')  =  {2«  +  l)[|J. 


(3.32) 


Figure  30  compares  the  zeroth  and  first  order  eigenvalues  to  the  results  obtained 
by  numerical  diagonalizaton  for  two  different  values  of  glco ,  0.05  and  0.25.  To 
facilitate  comparison,  the  values  actually  plotted  are  As  one  would  expect  from 

the  perturbative  calculation,  the  agreement  improves  for  larger  n  and  g,  as  Tio)(^l2 
becomes  less  significant  compared  to  the  other  terms  in  H.  In  addition,  note  how  the 
eigenvalues  naturally  form  two  sets  which  are  almost  smooth  functions  of  These 
two  sets  are  the  eigenvalues  corresponding  (in  the  appropriate  limit)  to  the  ±  signs  in 
Eq.  3.29. 


The  standard  result  for  first  order  correction  to  the  eigenvectors  is,  in  the 


notation  used  here. 


(3.33) 


Using  the  properties  of  the  displaced  number  states  and  a  substantial  amount  of 
algebra,  we  obtain  the  first  order  correction 


The  factorials  keep  the  distribution  narrow,  and  the  amplitudes  are  somewhat 

symmetrical  about  the  unperturbed  state.  To  use  these  states  in  calculation,  the  zeroth 

order  eigenstates  should  be  expressed  in  terms  of  the  atomic  eigenstates  and  photon 

number  states.  Although  not  conceptually  difficult,  the  calculation  is  cumbersome  and 

the  final  outcome  not  particularly  informative,  so  it  will  not  be  included.  However,  the 
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results  are  shown  in  Fig.  31,  where  the  probabilities  for  the  state  la,/?)  and  |/>,//)  (i.e.. 


the  atomic  energy  eigenstate  and  photon  number  state)  are  shown.  The  calculations 
for  this  figure  are  based  on  eigenstate  iT^i28+).  with  ^/ar=0. 17678. 

Let  us  now  return  to  the  goal  at  hand,  which  was  to  find  the  approximately 
factorizable  states  using  the  eigenstates  of  H.  Suppose  we  consider  the  outer  product 
of  an  atomic  state  given  by  the  polar  angle  6  in  the  Bloch  sphere  (let  the  azimuth  angle 
be  zero)  and  a  coherent  state  with  amplitude  a  =  -Jji ,  which  shall  be  denoted  \6,4^ ). 


In  general,  this  state  will  be  described  by  a  distribution  of  eigenstates  with  both 
positive  and  negative  parity.  The  parity  operator  of  Ref  26, 


O', 

a+— + 
0 


(3.35) 


is  easily  seen  to  have  as  eigenvectors  the  zeroth  order  eigenvectors  of  the  unperturbed 
Hamiltonian  (3.20).  Furthermore,  P  commutes  with  the  transformed  Hamiltonian 


H  =  +  T  -  is! <y)' 


(3.36) 


The  eigenvalues  of  like  parity  depend  approximately  linearly  upon  //  over  a 
restricted  range  of  n  (such  as  used  to  build  a  coherent  state)  as  shown  in  Fig  30.  The 
state  \0,^ )  is  composed  of  a  superposition  of  eigenstates,  which  will  have  time 


dependence 


=  (3.37) 

n 

± 


This  can  be  expressed  as 

I  2  I  Y,  I  )  .(3 . 3  8) 


Suppose  the  initial  state  has  been  chosen  so  that  the  C„.  coefficients  are  small  and  can 
be  neglected.  Now  by  expanding  as  a  Taylor  series  about  n=Jl ,  to  the  extent  that 
depends  linearly  on  //,  we  can  write 
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h'W)= 


ZC  exp 


Ti  3i 


(3.39) 


This  result  is  similar  in  form  to  the  Floquet  state  (3.2),  with  tio)'  = 


di 


Of 


course,  multiples  of  o)'  can  be  shifted  from  the  Floquet  exponent  to  within  the 
summation  without  effect  in  order  to  place  the  Floquet  exponent  in  the  desired  range. 

It  is  important  that  the  C„.  coefficients  can  be  neglected,  since  in  the  corresponding 
summation  for  the  minus  states  the  Floquet  exponent  and  <y' would  be  different. 
Compare  Fig.  32,  which  shows  the  variation  in  IXJco  using  the  eigenvalues  obtained 
from  the  first  order  perturbation  expression  given  by  Eq.  3.29,  to  Fig.  5.  In  the 
quantum  case  there  is  a  slight  dependence  upon  //,  but  otherwise  the  agreement  is  quite 
good.  Now  consider  the  inner  product  between  the  previously  mentioned  state 
1 6,  )  and  the  eigenstate  |  'fs- )  • 


(3.40) 


Finding  the  angle  ^  which  minimizes  A  will  produce  a  state  which  minimizes  the 
presence  of  the  minus  branch  (i.e.,  the  C„.’s)  and  therefore  remains  approximately 
factorizable.  Figure  33  shows  the  best  value  of  0for  various  values  of  e.  Two  sets  of 
eigenvectors  with  different  values  ofg/co  were  used  for  the  calculation,  0. 17678  and 
0.25.  The  agreement  with  Fig.  7  is  striking,  with  the  exception  of  small  e.  With  these 
values  of  glo),  the  small  s  region  requires  such  small  values  of  «  as  to  make  a 
quasiclassical  approximation  invalid.  The  figure  shows  the  agreement  improves  with 
decreasing  g/n>. 

Of  course,  some  presence  of  the  “wrong”  or  undesired  branch  is  unavoidable  in 
building  a  coherent  state.  However,  the  choice  of  0 given  above  reduces  it  to  a 
minimum. 
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We  now  have  found  two  different  techniques,  the  one  just  discussed  and  the  use 
of  Floquet  eigenstates,  which  lead  to  the  same  quasiclassical  states.  These  states  can 
be  combined  in  a  manner  which  leads  to  accurate  prediction  of  the  full  quantum  model. 
Let  us  now  use  these  results  with  the  understanding  of  chaos  in  the  semiclassical 
model  gained  in  Chap.  2  and  investigate  the  impact  of  semiclassical  chaos  on  the 
quantum  model. 
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CHAPTER  4 


CHAOS  AND  THE  QUANTUM  MODEL 

Finally  we  are  prepared  to  discuss  the  relationship  between  semiclassical  chaos 
and  the  behavior  of  the  fiilly  quantized  system.  Section  1 .3  introduced  some  of  the 
comparisons  between  semiclassical  and  quantum  models  done  by  other  authors.  In 
Sect  4. 1  we  will  focus  on  a  particular  approach,  examination  of  purity  and  uncertainty 
This  approach  has  been  used  extensively  by  the  authors  of  Ref  19-21  and  we  will 
show  some  of  the  pitfalls  which  may  trap  the  unwary.  In  addition,  some  problematic 
results  of  Graham  and  Hohnerbach  will  be  discussed.  At  last,  this  work  will  close  with 
Sec.  4.2,  in  which  the  effects  of  using  higher  values  of  angular  momentum  in  the 
atomic  system  are  considered. 

4.1  Impact  of  Semiclassical  Chaos 

Grigolini,  et  al.  have  published  a  series  of  papers’ ^'2'  on  investigations  of  the 
spin-boson  system,  both  with  and  without  the  RWA.  Their  approach  has  been  to 
compare  the  results  of  the  semiclassical  model  to  the  quantum  model  and  attempt  to 
find  correlations  between  the  onset  of  chaos  and  changes  in  dynamical  quantities, 
much  as  we  will  do  here.  In  studying  the  semiclassical  model,  they  generally  have  not 
focused  on  individual  trajectories,  but  performed  averages  over  an  initial  distribution 
which  satisfies  the  quantum  uncertainty  principle.  For  the  quantum  model,  they 
determine  the  evolution  using  the  Wigner  pseudoprobability  distribution  and 
Liouvillian  operators.  They  separate  the  Liouvillian  into  two  parts,  L^i  and  Lqqq.  L^i 
yields  the  same  results  as  the  semiclassical  model,  Eq.  2.4.  Lqq£3  (Quantization 
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Generating  Diffusion)  is  the  remainder  of  the  quantum  Liouville  operator  after 
removing  L(.|.  For  comparisons,  they  take  expectation  values  from  the  semiclassical 
distribution  and  from  the  corresponding  quantum  distribution.  This  is  a  rather  novel 
approach  and  holds  some  promise  for  understanding  the  relationship  between  the  two 
systems,  but  their  work  has  a  few  shortcomings  which  cast  doubt  on  their  conclusions. 
They  generally  consider  evolution  of  the  system  from  the  atom  in  the  lower  energy 
eigenstate,  and  do  not  fully  account  for  the  effects  of  the  two  quasiclassical  states 
(although  they  do  mention  one  effect,  the  splitting  of  the  field  along  the  plus  and  minus 
branches).  Furthermore,  they  use  the  RWA  results  as  a  basis  for  comparison,  and  any 
departures  from  the  RWA  results  they  attribute  to  semiclassical  chaos.  They  fail  to 
recognize  the  importance  of  initial  conditions  in  the  evolution  of  the  quantum  system, 
as  we  shall  see  below. 

An  example  will  serve  to  demonstrate  the  type  of  problem  which  can  arise  if  the 
presence  and  behavior  of  the  two  branches  is  not  recognized.  In  Ref.  20  they  suggest 
there  exists  a  qualitative  difference  in  the  dynamics  of  the  RWA  and  non-RWA 
quantum  models.  They  claim  this  difference  leads  to  a  dramatic  failure  of  the  RWA 
while  providing  good  agreement  between  the  quantum  and  semiclassical  non-RWA 
results.  Figure  3  of  Ref  20  plots  for  the  non-RWA  semiclassical  model  and  an 
approximation  they  call  the  Wigner  recipe,  while  in  Fig.  4  they  do  the  same  for  the 
RWA.  Parameters  as  defined  here  are  g==0.0898,  ojf=coq,  n=300,  and  f=0.778,  starting 
from  the  lower  energy  eigenstate  |A),  for  a  time  of  <y/=16;r.  The  Wigner  recipe 

consists  of  averaging  the  results  of  the  numerical  integration  of  the  semiclassical 
equations  (2.4)  over  an  initial  Wigner  distribution  for  a  coherent  state,  and  provides 
reasonably  good  agreement  with  the  full  quantum  model  for  the  calculations  shown. 

In  their  Fig.  4,  the  RWA  shows  a  collapse  occurring  about  t=S7r,  while  the  non-RWA 
quantum  calculation  shows  no  collapse  and  matches  the  semiclassical  model  very  well. 
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They  attribute  this  difference  to  a  failure  in  the  RWA.  It  is  more  correct,  however,  to 
attribute  it  to  comparing  two  different  initial  conditions.  For  the  given  values  of  the 
parameters,  |  ^+(0))  is  almost  coincident  with  |/)),  and  so  no  collapse  is  expected  and 

the  semiclassical  model  should  provide  good  agreement.  However,  in  the  RWA,  the 
state  \b)  is  formed  by  an  equally  weighted  linear  combination  of  the  plus  and  minus 

branches,  thus  the  collapse  is  expected.  A  more  appropriate  choice  would  be  to 
compare  the  results  along  the  plus  branch  in  both  instances,  in  which  case  the  RWA 
would  show  no  collapse  and  give  good  agreement  with  the  semiclassical  model. 
Alternatively,  one  could  compare  the  same  equally  weighted  linear  combination  of  plus 
and  minus  branches  in  the  non-RWA  as  used  in  RWA.  The  results  of  this  calculation 
are  shown  in  Fig.  34.  It  shows  a  collapse  in  time  of  about  /c=60,  compared  to  ^^=58. 5 
given  by  Eq.  3.13.  Thus,  the  main  difference  between  the  RWA  and  non-RWA  is  the 
length  of  a  small  difference  which  is  well  understood  from  the  material  in  Chap.  3. 

Now  let  us  turn  our  attention  to  Ref  21 .  In  this  article  they  develop  the 
Liouvillian  dynamics  and  separate  the  effects  of  classical  and  quantum  components.  In 
addition,  they  use  the  semiclassical  equations  of  motion  and  calculate  the  Lyapunov 
exponents  using  the  method  of  Ref  37,  as  was  done  here.  They  identify  a  threshold 
for  chaos  in  which  the  largest  Lyapunov  exponent  increases  to  almost  unity  from  about 
0. 1 .  They  refer  to  this  region  of  very  large  Lyapunov  exponent  as  the  chaotic  region, 
ignoring  the  region  with  smaller  Lyapunov  exponents.  This  is  not  typical  of  the  work 
in  this  field,  where,  for  example.  Ref  1 5  identifies  a  region  as  chaotic  for  a  Lyapunov 
exponent  as  small  as  0.087.  Likewise,  the  Henon-Heiles  system  was  considered 
chaotic  by  the  authors  of  Ref  37  with  a  Lyapunov  exponent  of  approximately  0.08. 

The  threshold  Grigolini,  et  al.  identified  for  highly  chaotic  behavior  occurs  at 
resonance  withg=0.5  and  17=5.  Because  of  the  small  photon  number,  their 
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conclusions  cannot  be  directly  applied  to  this  work,  since  such  a  small  photon  number 
precludes  the  quasiclassical  approximation  from  applying.  Nonetheless,  we  will 
continue  and  attempt  to  test  the  extension  of  their  work  into  the  parameter  space 
applicable  here. 

As  a  test  for  the  effect  of  semiclassical  chaos  on  the  quantum  system,  they  use 
the  entropy  of  the  atom.  Phoenix  and  Knight  present  a  method  for  finding  the  entropy 
using  Pa(,  the  atomic  reduced  density  operator. Given  the  matrix  elements 
and  pi,i„  the  eigenvalues  of  p^^  are 

«±  =  I {•  ±  -  Ai )'  +  4^61'  }  (41) 

and  the  entropy  is 

•^(Pat  )  =  -(«*  In  +  a.  In  a_ ).  (4.2) 

Using  the  expression  for  the  purity 

p=pl.+fL+2\pJ-  («) 

it  is  straightforward  to  relate  the  entropy  to  the  purity; 

.S’(P)  =  -^(l  +  V2P-l)ln[|(l  +  V2P-l)] 

-^(l-V2P-l)ln[i(l->/2P-l)] 

In  a  non-chaotic  region  of  the  parameter  space,  far  off-resonance,  Grigolini,  et 
al.  show  that  the  entropy  reaches  its  maximum  value  (In  2)  during  the  collapse  region, 
and  returns  briefly  to  zero  during  revivals.  They  then  claim  that  for  the  chaotic  system 
the  disappearance  of  the  revivals  in  the  quantum  model  is  a  result  of  the  chaos  in  the 
semiclassical  model.  The  absence  of  revivals  is  undeniable,  but  to  identify 
semiclassical  chaos  as  the  cause  is  rather  speculative,  as  we  will  show  momentarily. 

Let  us  turn  now  to  Ref  19,  which  reviews  the  material  in  the  previous  articles 
and  extends  it  somewhat.  Specifically,  they  add  a  discussion  of  the  entropy  growth 
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and  irreversibility  for  an  initial  photon  number  (Fock)  state  and  they  relate  the  growth 
of  field  uncertainty  to  semiclassical  chaos. 

The  entropy  growth  starting  from  a  number  state  is  included  in  order  to  compare 
with  results  from  the  RWA.  The  RWA  is  manifestly  non-chaotic  in  the  semiclassical 
model,  having  only  two  degrees  of  freedom.  In  spite  of  this  fact,  it  shows  an  overall 
increase  in  the  entropy  with  time  if  the  field  is  initially  in  a  coherent  state."'^  The 
increase  in  entropy  is  attributed  in  Ref  1 9  to  the  coherent  state  mimicking  a  heat  bath 
with  an  infinite  number  of  degrees  of  freedom.^®  The  collapses  result  from  the  many 
different  frequencies,  one  for  each  eigenstate  |<I>„)  of  the  oscillator.  Furthermore, 

because  the  frequencies  increase  as  ,  rather  than  integer  the  revivals  are  not 
perfect.  To  circumvent  the  irreversibility  and  consequent  increase  in  entropy  from 
starting  with  a  coherent  state,  in  Ref  1 9  they  calculate  the  entropy  with  the  field 
starting  in  a  number  state.  In  the  RWA,  if  the  field  begins  in  a  number  state  the  system 
evolves  reversibly,  while  in  the  non-RWA  it  does  not.  In  Ref  19,  Fig.  3  uses  the 
atomic  entropy  to  show  the  increasing  irreversibility  for  five  different  values  of  g,  with 
the  field  starting  in  a  number  state  and  the  atom  starting  in  the  lower  energy 
eigenstate.  They  claim  this  behavior  supports  their  conclusion  that  quantum 
irreversibility  is  a  signature  of  semiclassical  chaos.  However,  the  difference  between 
the  behavior  in  the  RWA  and  non-RWA  can  better  be  attributed  to  differences  in  the 
eigenstates  of  the  RWA  and  non-RWA  Hamiltonians.  The  eigenstates  of  the  RWA 
Hamiltonian  can  be  written 


(4.5) 


Thus  the  atom  in  the  upper  state  with  the  field  in  a  number  state  is  just 


purity 


As  shown  above  in  Sec.  3.3,  the  eigenstates  of  the  non-RWA  are  approximately 


non-RWA 
nt 


(4.7) 


The  forms  of  these  two  sets  of  states  are  very  similar.  Thus,  a  more  appropriate 
choice  for  the  initial  state  in  the  non-RWA  (in  order  to  compare  to  the  RWA)  would 
be 


|T(0)>  =  |+x)D'[^^j|/;),  (4.8) 

which  has  the  evolution  shown  in  Fig.  35.  The  calculation  shown  in  that  figure  should 
be  compared  to  Fig.  3e  of  Ref  19,  in  which  g=0. 898,  (resonance),  A7=10,  and 
£=1 .420,  and  which  shows  highly  irreversible  evolution.  In  contrast,  the  evolution  in 
Fig.  35  is,  to  a  good  approximation,  clearly  reversible,  even  when  the  semiclassical 
model  is  strongly  chaotic  and  the  approximation  leading  to  Eq.  4.7  is  not  particularly 
good.  This  comparison  again  demonstrates  the  importance  of  choosing  appropriate 
initial  conditions  when  comparing  the  results  from  the  RWA  and  non-RWA  models. 

Now  we  will  examine  the  quantum  uncertainty.  In  Ref  19  the  quantum 
uncertainty  is  defined  as 

(/  =  A(7, -Aa^.  (4.9) 

In  Fig.  4  of  that  reference,  U(J)  is  shown  for  five  values  of^{0.045,  0.090,  0.224, 
0.449,  0.898},  starting  with  the  field  in  a  coherent  state  with  /7=10  and  the  atom  in  the 
lower  eigenstate.  They  state  that  these  five  values  of  g  all  correspond  to  increasing 
values  of  the  Lyapunov  exponent  in  the  semiclassical  model,  and  that  figure  shows  that 
U{t)  increases  more  rapidly  for  larger  (hence  more  chaotic)  values  of  g.  However,  it  is 
difficult  to  attribute  the  growth  of  the  uncertainty  to  semiclassical  chaos  because  the 
smaller  values  of  g  show  growth  as  well  as  the  larger  values.  For  example,  their  Fig. 
4c-e  (g=0.224,  0.449,  0.898  respectively)  yield  the  same  general  behavior;  rapid 
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growth  in  U{t)  until  it  reaches  a  roughly  asymptotic  value,  which  is  larger  for  larger 
values  ofg^.  Curves  (c)  and  (d)  exhibit  the  greatest  similarity  to  each  other,  having 
roughly  the  same  slope  and  asymptotic  value.  Semiclassically,  however,  they  could 
not  be  more  different.  My  calculations  show  that  the  Lyapunov  exponent  for  the 
parameters  of  curve  (c)  is  converging  on  zero  and  is  clearly  non-chaotic  while  the 
parameters  of  curve  (d)  give  behavior  which  is  clearly  chaotic  with  a  Lyapunov 
exponent  of  approximately  0.09.  The  differences  in  the  evolution  of  these  two 
parameter  sets  semiclassically  and  similarities  quantum  mechanically  virtually  eliminate 
semiclassical  chaos  as  the  cause  of  the  uncertainty  growth. 

Grigolini,  et  al.  are  not  alone  in  mis-attributing  any  differences  observed 
between  the  RWA  and  non-RWA  behavior  to  chaos.  In  Ref  27  Graham  and 
Hohnerbach  do  the  same.  They  claim  the  disappearance  of  revivals  is  a  mark  of 
semiclassical  chaos,  as  are  changes  in  the  power  spectrum  of  the  field.  Figure  3  of 
Ref  27  shows  the  occupation  probability  of  the  upper  state  vs.  time  for  a  quantum 
trajectory  starting  from  the  state  \a)  withg/ry=0.01  and  17=100.  (The  occupation 

probability  of  the  upper  state  is  simply  (I+z)/2.)  The  trajectory  shows  definite  revivals 
and  collapses,  similar  to  Fig.  25.  They  contrast  that  to  the  behavior  shown  in  their  Fig. 
4,  for  which  glo}=0.075  and  all  else  is  the  same.  Their  Fig.  4  shows  no  clear  revivals, 
which  they  claim  is  the  result  of  semiclassical  chaos  based  on  their  analytic 
approximations.  Their  approximation  is  based  upon  studying  a  quantity,  K,  formed 
from  a  combination  of  spin  and  field  variables  which  is  conserved  in  the  RWA  but  not 
conserved  in  the  non-RWA.  Using  the  counter-rotating  terms  they  calculate  an 
expression  for  K,  which  they  evaluate  using  results  from  the  RWA.  This  process 
leads  to  a  prediction  that  the  entire  phase  space  will  be  chaotic  for  the  parameters  used 
in  their  Fig.  4.  Unfortunately,  these  parameters  (g/iy=0.075,  £=0.375,  77=100)  yield 
trajectories  which,  according  to  my  calculations,  are  completely  regular,  whether 
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starting  from  {a)  or  either  of  the  quasiclassical  states.  Furthermore,  the  absence  of 
revivals  is  not  surprising  because  the  given  value  of  £•  means  the  initial  state  |a)  is 
composed  mostly  of  the  minus  branch. 

They  find  an  analytical  result  for  the  power  spectrum  in  the  RWA  which  predicts 
two  peaks  symmetrically  displaced  from  nominal  field  frequency  o)  by  ±g/2^/jf .  These 
two  peaks  are,  of  course,  the  field  corresponding  to  the  plus  and  minus  branches. 

They  then  claim  the  asymmetry  due  to  the  dominance  of  one  peak  is  a  signature  of 
semiclassical  chaos.  Their  Fig.  6  is  used  to  support  this  claim,  where  they  plot  the 
power  spectrum  for  g/oj=0. 1 ,  £=0.4,  and  Ji =64.  This  plot  shows  two  peaks,  the  lower 
frequency  peak  larger  than  the  other  by  approximately  a  factor  of  2.9/0.35=8,  each 
centered  at  frequencies  ±0.002.  Using  the  Autler  Townes  model  we  find  the 
upper  atomic  state  is  composed  of  |a) «  0.30|  y/^)  +  0.95|  The  ratio  of  the  squares 

of  these  components  yields  approximately  9,  in  reasonable  agreement  to  the  ratio  of 
the  sizes  of  the  two  peaks.  Equation  2. 14  gives  an  interaction  induced  detuning  of 
0.002<y,  again  in  agreement  with  their  result.  And  finally,  like  the  previous  example 
from  their  paper,  the  semiclassical  model  is  not,  based  on  my  calculations,  chaotic  for 
this  set  of  parameters.  In  summary,  they  have  compared  their  results  of  quantum 
calculations  using  the  RWA  and  non-RWA,  attributing  differences  in  these  results  to 
the  presence  of  chaos  in  the  semiclassical  model.  The  differences  can  be  better 
explained  as  the  result  of  differences  in  the  quasiclassical  states  between  the  RWA  and 
non-RWA  models. 

What  then  can  be  said  about  the  correlations  between  semiclassical  chaos  and 
quantum  mechanical  evolution?  The  conclusion  of  the  present  work  is:  for  the  spin 
1/2  system,  little  or  nothing.  Recall  the  marked  transition  to  chaos  in  the  semiclassical 
model  for  the  plus  branch  near  Ji =90  with  £=1 .  Figure  36  compares  the  evolution  in 
the  Bloch  sphere  for  /» =100  and  /7=8 1 .  Semiclassically,  there  is  an  extensive 
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difference  between  the  trajectories.  The  chaotic  trajectory  covers  about  two-thirds  of 
the  surface  of  the  sphere,  while  the  regular  trajectory  remains  confined  to  the 
neighborhood  of  the  Autler-Townes  Floquet  eigenstate,  as  discussed  in  Sec.  2.3.  The 
relevance  of  the  quasiclassical  states  is  apparent,  since  the  quantum  trajectories  bear  a 
remarkable  resemblance  to  the  regular  semiclassical  trajectory.  In  addition,  the  two 
quantum  mechanical  trajectories  are  virtually  indistinguishable  from  each  other. 
Moreover,  the  similarity  between  the  two  trajectories  does  not  stop  at  the  Bloch 
sphere  evolution.  Figures  37-39  compare  the  purity,  spectra,  and  uncertainty  product 
for  the  two  quantum  trajectories. 

In  Fig.  37,  one  main  difference  in  the  purity  for  the  two  quantum  trajectories  can 
be  readily  observed.  For  the  trace  for  /7=8 1 ,  a  few  oscillations  are  visible  at  the  very 
beginning  and  again  near  /=550.  These  correspond  to  a  partial  collapse  and  revival 
due  to  the  unavoidable  presence  of  some  portion  of  the  minus  branch  field,  as 
discussed  in  Sec.  3.3.  The  partial  collapse  and  revival  are  also  visible  in  the  trace  for 
7=100,  but  around  /=675.  The  ratio  of  these  two  times,  according  to  Eq.  3.14, 
should  be  equal  to  the  ratio  of  7  for  each  trajectory.  In  fact,  the  ratio  is  approximately 
0.8 1 4,  in  agreement  with  the  predicted  value  of  0.8 1 .  However,  there  are  no  other 
significant  differences  between  the  trajectories  to  indicate  the  semiclassical  system  has 
undergone  a  transition  to  chaos. 

A  similar  comparison  can  be  made  from  Fig.  38.  The  spectra  from  the  quantum 
trajectories  are  shown  together  in  Fig.  38a.  The  major  features  of  each  spectrum  are 
identical.  The  differences  in  the  two  spectra  at  the  base  of  some  of  the  peaks  are 
artifacts  of  the  Fourier  transform  resulting  from  the  difference  in  the  interaction 
induced  detuning  for  the  two  trajectories.  On  the  other  hand.  Fig.  38b  contrasts  the 
two  spectra  for  the  semiclassical  system.  For  7=100  two  small  peaks  can  be  seen  near 
the  odd  values  of  <y,  while  for  7=8 1  these  are  regions  of  noise.  Likewise,  the  base  of 
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each  even  value  of  (o  is  a  region  of  noise  for  the  chaotic  trajectory.  The  two  small 
peaks  in  the  «  =  I00  spectrum  result  from  the  Autler-Townes  model  not  giving  the 
precise  starting  point  needed  for  a  periodic  trajectory  in  the  semiclassical  model.  By 
trial  and  error,  though,  one  can  “zero  in”  on  the  optimum  starting  point  and  make  the 
pair  as  small  as  desired.  However,  trial  and  error,  even  to  8  decimal  place  accuracy, 
does  not  get  rid  of  the  noise  in  the  chaotic  spectrum  or  provide  a  stable  trajectory. 

Finally,  Fig.  39  shows  the  evolution  of  the  quantum  uncertainty,  f/,  given  by  Eq. 
4.9.  In  this  figure,  (/is  plotted  against  a  scaled  time,  tin,  in  order  to  suppress  the 
difference  in  the  two  trajectories  due  to  the  interaction  induced  detuning.  Without  the 
scaling  the  uncertainty  for  «  =81  grows  somewhat  faster  than  for  «=100.  However, 
the  scaling  demonstrates  that  the  faster  growth  is  not  due  to  chaos,  but  to  the 
interaction  induced  detuning  as  the  field  phase  uncertainty  grows  (see  Fig.  26).  The 
two  trajectories  do  not  significantly  differ  until  about  one  half  of  the  revival  time,  when 
the  coherent  state  approximation  begins  to  fail.  Even  then,  the  difference  in  the  two 
trajectories  is  small. 

To  further  illustrate  the  absence  of  a  link  between  semiclassical  chaos  and  the 
behavior  of  the  quantum  model,  consider  now  the  evolution  along  the  minus  branch. 
The  next  series  of  figures.  Figs.  40-43,  starts  along  the  minus  branch,  which  does  not 
display  chaos  for  «  >30.  The  main  differences  between  these  trajectories  and  those  of 
the  previous  set  of  figures  are  in  the  semiclassical  trajectory  for  «  =81.  Figure  40d 
shows  that  the  orbit  for  n  =8 1  remains  stable,  while  Fig  42b  shows  that  the  spectra  for 
w=81  and  /7=100  are  virtually  identical.  The  quantum  model  trajectories  differ  only 
slightly  along  the  plus  and  minus  branches.  The  minus  branch  trajectories  are  less  pure 
than  their  plus  branch  counterparts,  as  shown  in  Fig.  41.  In  addition,  the  revivals  are 
even  more  obscured  for  the  minus  branch.  These  observations  run  directly  counter  to 
the  contentions  of  Grigolini  and  Graham,  since  the  minus  branch  is  not  chaotic. 
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Likewise,  Fig.  43  shows  that  the  uncertainty  product  for  the  minus  branch  trajectories 
is  higher  than  the  uncertainty  product  for  the  plus  branch  trajectories.  Although  it  is 
notable  that  the  uncertainty  products  are  closer  in  value  for  the  two  minus  branch 
trajectories,  the  larger  value  for  the  minus  branch  shows  there  is  no  correlation  with 
the  chaotic  behavior  of  the  semiclassical  model.  These  examples  demonstrate  the 
absence  of  a  direct  correlation  between  semiclassical  chaos  and  the  behavior — whether 
measured  by  purity,  entropy,  spectrum,  or  uncertainty — of  the  quantum  model. 
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4.2  Higher  values  of  Angular  Momentum 


The  difficulty  in  establishing  a  correspondence  between  the  onset  of  chaos  in  the 
semiclassical  problem  and  a  signature  in  the  quantum  problem  seems  to  be  related  to 
the  highly  non-classical  nature  of  a  two  level  atom.  The  field  is  as  classical  as  possible 
in  a  coherent  state,  which  has  the  minimum  possible  uncertainty.  However,  the  spin 
state  can  be  made  more  classical  by  moving  to  higher  values  of  angular  momentum 
quantum  number,/  Increasingy  corresponds  to  an  interaction  between  the  field  and  a 
larger  number  of  atoms,  N,  where  N=2j.  Graham  and  Hdhnerbach  showed  the 


distribution  of  eigenvalues  becomes  typical  of  a  chaotic  system  for  larger  values  of 
y  25,27  Moving  to  higher  values  of  j  might  make  the  quasiclassical  states  evolve  more 
closely  to  the  semiclassical  states  in  the  chaotic  regime. 

Numerically,  calculating  trajectories  for  higher  values  of  j  is  much  more 
challenging:  The  larger  number  of  atomic  states  to  calculate  is  compounded  by  the 


scaling  given  immediately  following  Eq. 


n—>n2J 


^=MlJ 

y=M/j 


2.8; 


(4.10) 


principally,  that  n  increases  linearly  with  j.  Increasing  from  y=l/2  to  y=3/2  results  in  a 
roughly  fivefold  increase  in  computations  at  each  time  step.  Furthermore,  increasing 
n  results  in  higher  frequencies  which  must  be  integrated  using  smaller  time  steps.  For 
the  first  example  given  below,  increasing  from y=l/2  to 7=3/2  resulted  in  an  increase  in 
computation  time  by  a  factor  of  forty.  For  this  reason,  this  investigation  was  halted  at 
7=3/2,  even  though  higher  values  of  j  might  prove  to  be  interesting. 
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The  Hamiltonian  (1 . 1)  for  the  higher  j  states  becomes 
H  =  hcOf)J^  +  tioxi^a  +  Tig{j^  +  J.  ){a  +  a^). 


(4.11) 


(4.12a) 


(4.12b) 


(4.12c) 


Obviously,  the  eigenstates  of  are  the  basis  states  for  describing  the  state  of  the  three 
atoms.  While  the  state  of  the y=l/2  system  can  be  completely  specified  by  two 
parameters,  such  as  polar  angles  6  and  (j)  (up  to  an  unimportant  overall  phase),  this  is 
not  the  case  for y=3/2.  Instead,  the  equivalent  initial  condition  is  specified  as  an 
eigenstate  of  the  operator  J  -n,  where  ti  is  a  unit  vector  specified  by  polar  angles  6 
and  (f).  Using  the  standard  angular  momentum  relationships  J=iJj+J^I2  and 


)/2/,  (quite)  a  bit  of  algebra  yields  the  starting  points  for  the  plus  and  minus  branches 
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These  are  the  eigenstates  of  J  -n  with  eigenvalue  ±3/2,  and  become  the  special  plus 
and  minus  branch  starting  points  with  0=0^. 


There  are  no  tremendous  differences  in  the  evolution  of  the  quasiclassical  states 
for  these  values  of y,  although  some  interesting  effects  do  begin  to  arise.  Figure  44 
shows  the  Bloch  sphere  evolution  for  the  small  ^case  shown  in  Fig.  18,  but  for  a 
shorter  time  (/=126).  The  parameters  used  in  the y=3/2  model  are  £=0. 141,  77=96,  and 
^0.02887  (0.05/V3),  which  correspond  to  £=0. 141,  77=32,  and  ^0.05  for y=l/2. 

The  only  significant  difference  is  that  the  deviations  from  the  Autler-Townes  model 
(which  are  just  small  amplitude  Rabi  oscillations)  at  the  beginning  of  the  trajectory  are 
smaller.  In  this  way  the  larger  j  model  does  match  the  semiclassical  model  more 
closely  than  the  j=\l2  model.  Since  both  of  these  sets  of  parameters  scale  to  the  same 
parameters  for  the  semiclassical  model,  we  expect  the  detuning  will  be  the  same  for 
both  values  of y,  which  is  indeed  the  case,  based  on  the  complex  phasor  method  of 
finding  the  detuning. 

As  fis  increased,  the  differences  between  the y=l/2  and  the y=3/2  models 
become  somewhat  more  pronounced,  and  the y=3/2  model  begins  to  show 
characteristics  which  are  reminiscent  of  the  semiclassical  model.  Let  us  examine  these 
characteristics  using  a  set  of  trajectories  for  the y=3/2  model  which  have  parameters 
corresponding  to  those  already  seen,  Sr=\  and  77=81  or  77=100.  For  the y=3/2  model 
the  corresponding  values  are  s=\  and  77=243  or  77=300  The  differences  between  the 
y=l/2  and  the y=3/2  models  are  first  observable  in  Fig.  45,  in  which  we  can  see  the 
trajectories  begin  to  spiral  inward  as  the  system  loses  purity  more  quickly  than  for  the 
y=I/2  model  (refer  to  Fig.  36  and  Fig.  40),  especially  in  regard  to  the  plus  branch.  The 
difference  between  the  plus  and  minus  branches,  which  was  almost  indiscernible  for 
y=l/2,  is  now  quite  distinct. 

Another  difference,  although  subtle,  is  demonstrated  by  Fig.  46.  This  figure 
shows  the  Bloch  sphere  from  a  different  orientation,  looking  along  the  +y  axis.  In 
three  of  the  plots  (a),  (c),  and  (d),  the  trajectories  are  roughly  symmetrical  about  x=0. 
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However,  in  (b)  the  trajectory  is  shifted  slightly  in  the  +x  direction.  This  particular 
plot  corresponds  to  the  chaotic  semiclassical  trajectory  shown  in  Fig.  36d,  which 
clearly  shows  a  strong  shift  toward  +x.  Figure  47  contrasts  x{t)  for  the  plus  branch  in 
the 7=3/2  and y=l/2  models.  Both  trajectories  correspond  to  £=1.0  and  /7=81  in  the 
semiclassical  model.  In  addition  to  the  faster  oscillations  with  ry=(2«+l)<y',  x{i)  for  the 
7=1/2  model  exhibits  normal,  slower  Rabi  oscillations  under  a  Gaussian  collapse 
envelope.  These  are  similar  to  the  oscillations  observed  in  the  Autler-Townes  model 
given  by  Eq.  2.9  for  Y;!i0.  Of  course,  in  the  Autler-Townes  model  the  oscillations  do 
not  collapse.  However,  the  oscillations  for  the 7=3/2  model  are  different.  Figure 
47shows  the  maximum  value  of  x  reached  on  each  fast  oscillation  changes  very  little, 
while  the  minimum  value  of  x  increases  considerably.  This  behavior  results  in  the  shift 
toward  +x,  not  unlike  the  shift  seen  in  the  semiclassical  model,  and  is  perhaps  a 
signature  of  the  transistion  to  chaos. 

The  purity  also  shows  a  possible  signature  of  the  semiclassical  chaos,  as  shown 
in  Fig.  48.  While  for 7=1/2  the  purity  of  the  (non-chaotic)  minus  branch  was  lower 
than  the  purity  for  the  (chaotic)  plus  branch,  for 7=3/2  the  purity  of  the  plus  branch  is 
lower.  It  is  not  unexpected  that  the  purity  along  a  particular  branch  will  be  slightly 
lower  for  smaller  JJ,  as  observed  for  the  minus  branch.  However,  the  purity  for  the 
plus  branch  is  substantially  less  for  the  smaller  JJ  trajectory,  again  suggesting  the  lower 
purity  could  conceivably  be  a  signature  of  semiclassical  chaos. 

Finally,  the  uncertainty  product  for  the  j=3/2  model  yields  results  similar  to  those 
for  the  purity;  The  uncertainty  product  is  significantly  larger  for  the  plus  branch  than 
for  the  minus  branch,  a  reverse  of  the 7=1/2  result.  However,  Fig.  49  shows  the 
uncertainty  product  is  almost  independent  of  JJ  along  a  particular  branch  when  plotted 
against  the  scaled  time. 
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Clearly  more  work  can  be  done  in  this  area.  The  results  here  are  intriguing 
because  they  suggest  the  possibility  the  signatures  of  chaos  put  forward  by  Grigolini, 
et  al.  may  be  applicable  for  higher  values  of  j.  However,  these  results  do  not  support 
a  definitive  statement  for  either  case.  Additional  work  must  be  done  by  extending  the 
integration  time,  examining  the  transition  to  chaos  at  other  points  in  the  parameter 
space,  studying  the  eigenvalue  structure,  and,  finally,  repeating  the  process  for  higher 
values  of  j.  Although  the  resources  were  not  available  to  undertake  these  additional 
investigations  in  support  of  the  present  effort,  the  initial  investigation  has  opened  the 
possibility  for  further  work  in  this  field. 


114 


Figure  44.  Bloch  sphere  evolution  for  the 7=3/2  model  with  £=0. 141  and  n=96. 
Parameters  correspond  to  the  semiclassical  parameters  £=0. 141  and  Ti  =32,  the  same  as 
_ Fig.  18,  for  a  time  of /=50.  (a)  Schrodinger  picture,  (b)  Interaction  picture. 
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4.3  Conclusion 


This  work  began  with  a  review  of  the  recent  and  not-so-recent  efforts  to 
investigate  the  behavior  of  the  spin-boson  system.  Much  work  has  been  done  using 
various  approximations  with  the  field  treated  both  classically  and  quantum 
mechanically.  With  a  semiclassical  field  the  system  had  been  shown  to  exhibit  chaotic 
behavior  in  the  realm  of  large  numbers  of  atoms  and  an  initially  weak  field.  We  then 
studied  some  of  the  tests  for  chaos  and  found  another  region  of  chaos:  single  atoms 
and  intense  fields.  A  special  set  of  states  which  evolved  quantum  mechanically  and 
semiclassically  in  close  agreement  were  identified  and  studied.  In  spite  of  their 
similarity,  one  of  these  trajectories  became  unstable  before  the  other.  The  reason  was 
identified  as  a  merging  of  harmonics  in  the  spectra  which  resulted  in  the  trajectories 
becoming  unstable.  Meanwhile,  the  quantum  trajectories  did  not  change  their 
behavior.  Thus  we  have  seen  clearly  an  example  in  which  the  semiclassical  system 
undergoes  a  transition  from  regular  to  chaotic  behavior,  while  the  quantum  system  is 
qualitatively  unchanged.  This  example  strongly  supports  the  contention  that  the 
changes  observed  in  the  purity  (or  entropy)  and  the  quantum  uncertainty  are  unrelated 
to  the  presence  of  chaos  in  the  semiclassical  model.  By  carefully  comparing  the 
semiclassical  and  quantum  models  using  the  quasiclassical  states,  we  have  found  there 
is  no  essential  difference  in  the  quantum  evolution  between  the  manifestly  non-chaotic 
RWA  and  the  non-RWA  models,  so  long  as  the  quasiclassical  approximation  holds. 
The  quasiclassical  approximation  holds  until  ejn  (or  gj -Jh  )  becomes  too  large  or  77 

becomes  too  small.  Most  importantly,  the  key  in  comparing  the  RWA  to  the  non- 
RWA  is  in  comparing  similar  states — whether  comparing  evolution  along  a  particular 

branch  or  along  a  linear  combination  of  branches. 
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Finally,  an  intriging  area  has  been  opened,  that  of  the  dynamics  of  the  system  for 
larger  values  of  angular  momentum  quantum  number,  y.  The  one  example  examined 
here,y=3/2,  indicated  the  possibility  that  certain  signatures  of  chaos  may  be  present. 
The  quasiclassical  trajectories  for  this  value  of  j  were  not  more  factorizable  in  the 
region  near  chaos.  However,  it  did  show  share  some  characteristics  with  the 
semiciassical  model,  and,  as  in  all  fields  of  research,  more  work  remains  to  be  done. 
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APPENDIX  A 


Matrices  used  in  Perturbation 
Calculation 

In  Sec.  2.3  a  calculation  to  find  the  frequencies  which 
appear  in  the  semiclassical  model  was  presented.  However, 
the  matrices  were  sufficiently  large  to  be  impractical  for 
inclusion  in  the  main  text.  Therefore  they  have  been 
relegated  to  this  appendix.  Equation  2.20  is  a  sum  of  matrix 
equations,  but  because  some  of  the  matrix  elements  contain 
complex  exponentials,  the  terms  in  the  sum  are  not 
independent.  Therefore,  the  summation  needs  to  be 
rewritten  into  a  single  matrix  equation  which  can  be 
numerically  and  perturbatively  solved.  The  basis  vectors 
can  be  broken  into  five  parts  corresponding  to  Sa^,  Sa2,  5x, 
Sy,  and  &.  Each  of  these  five  parts  has  an  infinite  number 
of  elements  corresponding  to  the  coefficients  of  the 
complex  exponential  An  arbitrary  vector  x  may  be 

written  as  shown  to  the  right.  Furthermore,  suppose  we 
designate  the  subvectors  of  x  corresponding  to  the  different 
expectation  values  as  x^j,  x^  x^,  etc.  Likewise,  we  can 
designate  submatrices  of  the  matrices  A  and  B  as  A^^, 

Ba,,7,  etc.,  using  the  standard  row,  column  index 
convention.  Now  we  can  determine  specifically  the  various 
non-zero  parts  of  A.  The  simplest  parts  are  those  identified 
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by  ±1  in  Eq.  2.20,  These  submatrices  are  merely  the 

identity  matrix,  with  ±1  along  the  main  diagonal  of  the  submatrix.  The  next  simplest 


to  determine  are  the  submatrices  along  the  main  diagonal  of  A.  These  all  have  the 
form 
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^X.x  ^y.y  ^2.2 


Now  comes  the  cosine  expansion.  Because  of  the  two  exponentials  forming  the 
cosine  function,  the  index  of  each  element  multiplied  by  this  term  will  be  raised  or 
lowered  by  1,  thus 
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=  -A 


The  first  part  of  the  matrix  B  is  also  straightforward.  The  submatrix  ^  is  the 
identity  matrix  multiplied  by  the  constant  -ji.  To  obtain  the  next  part  we  will  return  to 
the  equation  from  which  the  matrices  are  extracted,  considering  only  the>'  portion. 
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=  —z 


(A.4) 


Now  expanding  as  a  sum  of  complex  exponentials 


n,m=-oo 

M~n  i,n 


(A.  5) 


From  this  expression  it  is  evident  that  the  matrix  element  ^  ^  =  -z^^„,  or 


(A.6) 


Likewise,  B^^^  will  have  the  same  form,  except  with  in  place  of  z. 

In  practice,  the  submatrices  were  each  truncated  at  3 1  rows  to  calculate  the 
unperturbed  eigenvectors,  for  a  total  of  155  elements  in  each  eigenvector  and  24,025 
elements  in  each  matrix. 
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APPENDIX  B 


A  Note  on  Numerical  Methods 

A  large  portion  of  this  research  has  been  done  numerically,  and  so  a  note  on  the 
methods  employed  is  in  order.  Several  types  of  calculations  were  accomplished 
numerically:  algebraic  calculation,  integration,  and  eigensystem  solution.  Most  of  the 
computations  were  performed  using  a  personal  computer  (40  MHz  80386SX  with 
80387  coprocessor  and  8  MB  memory),  although  some  were  performed  on  the 
University  of  Arkansas  mainframe  (Sun  SparcCenter  2000  with  eight  40  MHz 
microprocessors  and  256  MB  memory). 

Algebraic  calculations  were  generally  accomplished  using  Mathcad.  This 
includes  evaluations  of  the  Floquet  exponent,  the  X’s  and  B's,  etc.  Analysis  of  the 
numerical  integrations  (plots  of  time  series,  FFTs,  Poincare  sections,  etc.)  was  also 
done  in  Mathcad.  Most  plots  in  the  main  text  are  Mathcad  output.  In  addition,  some 
of  the  symbolic  manipulation  capabilities  of  Mathcad  were  used  during  analytical 
work. 

The  numerical  integration  was  accomplished  using  a  fourth-order  Runge-Kutta 
routine  with  adaptive  step  size  control.  The  algorithm  and  much  of  the  code  was 
taken  from  the  book  Numerical  Recipes  in  Pascal  (Press,  Flanneiy,  Teukolsky,  and 
Vetterling,  Cambridge  University  Press,  1989).  The  adaptive  step  size  control  is 
implemented  by  taking  a  single  step  in  time,  then  repeating  the  calculation  in  two 
steps,  each  using  one  half  the  original  step  size.  If  the  difference  in  the  result  is 
sufficiently  small,  then  the  result  is  saved  and  a  new  step  size  is  computed  based  on  the 
difference.  If  the  difference  is  too  large,  a  new  step  size  is  computed  and  the  process 
repeated.  The  accuracy  parameter  {eps  in  the  reference  given  above)  was  set  at  10'^ 
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to  achieve  a  balance  of  execution  speed  and  accuracy.  The  program  on  the  PC  was 
implemented  using  Turbo  Pascal  for  Windows.  Semiclassical  calculations  required 
about  2  minutes  for  C=200  in  the  examples  used  in  Chapter  4.  Quantum y=l/2 
calculations  required  about  2  hours  for  the  corresponding  trajectories.  All  the 
numerical  integration  for  the y=3/2  model  was  done  on  the  mainframe.  The  same 
integration  routine  was  used.  Using  the  parameters  of  the  first  example  in  Sect.  4.2, 
for  /=4;t,  the y=l/2  case  required  about  2  minutes  of  real  time  to  complete,  while  the 
7=3/2  case  required  about  80  minutes  of  real-time.  CPU  time  was  approximately  one 
half  of  those  amounts.  The  later  examples  for  the 7=3/2  case  required  15-16  hours  of 
CPU  time  to  integrate  up  to  (=50.  The  normalization  of  the  waveflmction  and  the 
total  energy  were  computed  to  check  the  accuracy  of  the  integration,  and  both  were 
typically  conserved  to  1  part  in  10^-10^. 

The  eigensystems  were  solved  using  a  combination  of  programming  and 
Mathcad.  A  simple  program  from  Numerical  Recipes  in  Pascal  was  employed  to 
calculate  the  eigenvalues  of  (1 . 1),  which  is  tridiagonal  if  separated  into  its  even  and 
odd  parity  parts.  All  other  eigenvalues  and  all  eigenvectors  were  computed  using  the 
built-in  eigenvalue  and  eigenvector  functions  of  Mathcad.  The  calculation  of 
eigenvalues  typically  required  a  few  minutes  to  complete  on  the  PC,  while  the 
calculation  of  a  full  set  of  eigenvectors  could  take  several  hours. 
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